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HIGHER HARMONIC CONTROL ANALYSIS FOR VIBRATION 
REDUCTION OF HELICOPTER ROTOR SYSTEMS 

Khanh Q. Nguyen 
Ames Research Center 


SUMMARY 

An advanced higher harmonic control (HHC) analysis has been developed and applied to investigate 
its effect on vibration reduction levels, blade and control system fatigue loads, rotor performance, and 
power requirements of servo-actuators. The analysis is based on a finite element method in space and 
time. A nonlinear time domain unsteady aerodynamic model, based on the indicial response formulation, 
is used to calculate the airloads. The rotor induced inflow is computed using a free wake model. The 
vehicle trim controls and blade steady responses are solved as one coupled solution using a modified 
Newton method. A linear frequency-domain quasi-steady transfer matrix is used to relate the harmonics 
of the vibratory hub loads to the harmonics of the HHC inputs. Optimal HHC is calculated from the 
minimization of the vibratory hub loads expressed in term of a quadratic performance index. Predicted 
vibratory hub shears are correlated with wind tunnel data. The fixed-gain HHC controller suppresses 
completely the vibratory hub shears for most of steady or quasi-steady flight conditions. HHC actuator 
amplitudes and power increase significantly at high forward speeds (above 100 knots). Due to the applied 
HHC, the blade torsional stresses and control loads are increased substantially. For flight conditions 
where the blades are stalled considerably, the HHC input-output model is quite nonlinear. For such 
cases, the adaptive-gain controller is effective in suppressing vibratory hub loads, even though HHC 
may actually increase stall areas on the rotor disk. The fixed-gain controller performs poorly for such 
flight conditions. Comparison study of different rotor systems indicates that a soft-inplane hingeless 
rotor requires less actuator power at high speeds (above 130 knots) than an articulated rotor, and a 
stiff-inplane hingeless rotor generally requires more actuator power than an articulated or a soft-inplane 
hingeless rotor. Parametric studies for a hingeless rotor operating in a transition flight regime and for an 
articulated rotor operating at the level-flight boundary (high speed and high thrust conditions) indicate 
that blade parameters including flap, lag, torsion stiffness distributions, linear pretwist, chordwise offset 
of center-of-mass from elastic axis and chordwise offset of elastic axis from aerodynamic center can be 
selected to minimize the actuator power requirements for HHC. 



1 INTRODUCTION 


Inherent in the aircraft’s design, vibration has remained one of the major problems facing helicopter 
development for years. With current passive vibration alleviation methodologies, vibration in helicopters 
has been reduced to an acceptable level over most of the flight envelope at the price of a substantial 
weight penalty. Currently, the maximum speed and maneuvering capabilities for most of the modern 
helicopters are limited by excessive vibration. With growing demands to expand the helicopter flight 
envelope and stringent requirements for lower vibration levels, effective and efficient vibration control 
systems are required. In this regard, the higher harmonic control systems have the potential for giving 
the helicopter a “jet-smooth” flight. 

1.1 Helicopter Vibrations 

Unlike fixed-wing aircraft, where vibrations primarily originate from the engines or are caused by 
atmospheric turbulence, the major source of vibration in a helicopter is the rotor. In fact, helicopter 
vibration is an aeroelastic problem which involves the complex interactions between the aerodynamic, 
inertial, and structural loads. The helicopter rotor, designed primarily to hover, is also required to have 
a forward flight capability and so must be able to move edgewise through the air. The primary source 
of helicopter vibration is the aerodynamic environment on the rotor disk. In this environment, a rotor 
blade section encounters a periodic variation in the free-stream velocity. At high forward speeds, the 
blades experience substantial transonic effects at the advancing tips, dynamic stall and reverse flow on 
the retreating side, and swept flow in the fore and aft regions of the rotor disk. Due to the aerodynamic 
environment and the inherent blade flexibilities, the rotor blades undergo large periodic motions which, 
in turn, are fed back as a source of input driving the aerodynamic loads. Therefore, rotor analysis must 
be considered as an aeroelastic problem since the blade motions are coupled with the aerodynamic loads. 

In addition, the temporal and spanwise variations in the blade airload forms in a system of shed and 
trailed vorticities behind the rotor disk. This vortex system induces an unsteady rotor nonuniform inflow 
which directly affects the blade angle of attack and, in turn, influences the blade airloads. For certain 
flight conditions, close blade-vortex interaction is also a significant source of helicopter vibration. As in 
the case of the blade motion-airload interaction, the induced inflow-airload interaction forms a second 
closed-loop in the rotor aeroelastic problem. 

The periodic variations in blade free-stream velocity and blade angle of attack produce blade 
airloads at many harmonics of the rotor rotational speed. Combined with the higher harmonic inertial 
loads, these result in multiharmonic blade loads transmitted through the rotor hub to the airframe and 
hence as a source of helicopter vibration. 

In general, not all harmonics of the blade loads are transmitted to the fuselage. This is certainly 
true for a properly tracked rotor with equally spaced blades. For a tracked rotor, all blades have identical 
aeroelastic behavior except for an equal phase difference around the rotor azimuth. The blade loads 
consist of many harmonics of the rotor rotational speed which tend to lower magnitudes at the higher 
frequencies. For the vertical hub force and the rotor torque, only the blade loads harmonics at pN b per 
rev (where N b is the blade number, p is a nonnegative integer) are summed to yield the hub loads at 
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pN b per rev. For the inplane hub shears and moments, only the blade load harmonics at pN b - 1 and 
pN b + 1 per rev contribute to the fuselage vibration at pN b per rev. The sum of all other harmonics 
cancels out at the hub. 

Regardless of which harmonics of the blade loads are transmitted through the hub, the resultant 
hub loads in the nonrotating frame — or fixed frame — occur only at pN b per rev. For the case when 
V is zero, the steady components of the hub loads act as a source of steady forces and moments on 
the fuselage, while the hub loads harmonics at N b per rev generally contribute to airframe vibration. 
The higher harmonics of hub load at 2N b ,3N b , . . . per rev are typically much smaller than the N b 
per rev components and thus can be neglected. Although these vibratory hub loads dominate the 
vibration problem, they are not the only source of helicopter vibration. The vibration contributions from 
the engine, wake induced loads acting on the fuselage and control surfaces, the tail rotor, and small 
differences between blades are typically small. 

Vibration has remained one of the major problem for helicopters. It affects the helicopter handling 
qualities, contributes to the fatigue of structural components, reduces the reliability of on-board electronic 
equipment, and influences weapon aiming precision in military aircraft. Vibration affects ride comfort, 
adds to the fatigue of pilots, crew, and passengers and also increases maintenance time and cost. The 
high vibration levels experienced by a helicopter could in many cases pose a limitation to the vehicle 
forward speed and maneuvering capabilities. 

Currently, there are several approaches to reduce helicopter vibration (refs. 1-3). A natural ap- 
proach, usually performed at the vehicle design phase, involves the placements of blade and fuselage 
natural vibration frequencies away from the rotor harmonic frequencies. These tuning effects help avoid 
blade and fuselage resonances with the rotor harmonics. This method, however, requires extensive mod- 
eling and testing of the blade and fuselage structural characteristics throughout the aircraft development 
phase. Recently, structural optimization techniques are under consideration to design the rotor and the 
fuselage in order to reduce the transmitted vibratory hub loads and the fuselage vibration. 

In spite of careful rotor design, significant vibrations exist in flight, and this necessitates the addition 
of further vibration control devices. The passive vibration control devices currently used in helicopters 
can be classified according to whether the vibrations are reduced at the rotor, at the fuselage, or by 
isolating the rotor from the fuselage. Dynamic absorbers are the vibration control devices mounted 
on blades or rotor hub. These devices can be considered essentially as spring-mass systems and are 
tuned to vibrate in such a way that the oscillatory blade loads are canceled before being transmitted to 
the rotor hub. The problem of pendulum absorbers is that each mode of vibration requires a separate 
absorber, and this makes the complete absorber system rather complicated. Also, due to the inability of 
current rotor analysis to accurately predict the vibratory loads, pendulum absorbers can only be tuned 
through repetitive flight tests. 

Rotor isolation is another popular vibration reduction technique and has been applied in several 
variations. Flexible mounting of the rotor structures to the fuselage can be designed in such a way as 
to minimize the transmission of the vibratory hub loads. In another variation, the rotor and the engine 
gearbox are mounted on a nodal-beam structure. The beam structure is configured such that the vehicle 
airframe is connected to the beam nodal points; hence, under excitation, only small oscillatory forces 
are transmitted to the fuselage. Another isolation device, typically known as the dynamic antiresonant 


3 



vibration isolator (DAVI), combines some of the best features of isolation and dynamic absorbers. 
Basically, this device consists of spring elements to which pendulums are attached. With proper tuning, 
the work done by the vibratory hub loads will be absorbed by the spring mechanism strain energy 
or consumed by the pendulum kinetic energy. In general, DAVI offers much better performance than 
simple isolation or absorber devices. 


The vibration reduction methods in the fuselage use either vibration absorbers which are tuned to 
the troublesome frequencies and placed at selected locations in the fuselage or vibration isolators in 
the form of soft mounting of the cabin floor to the airframe structures. With proper airframe design, 
the critical locations can be placed close to the nodes of the fuselage excited modes. In general, the 
methods of reducing vibration in the fuselage can only achieve local effects. 

Together, these four techniques represent the vibration control technology currently used on pro- 
duction helicopters. With proper design, tuning, and selection, they can provide a relatively comfortable 
ride over a large part of the helicopter flight envelope. However, these devices contribute to a significant 
weight penalty for the aircraft. Furthermore, their performance is quite restrictive and can deteriorate 
at high speeds and during maneuvers. Also, they are typically tuned only to selected locations in the 
airframe and work effectively for a certain range of the aircraft gross weight, payload distribution, and 
fuel loading. 

Up to now, the specifications for the desired vibration levels have generally been lower than the 
levels that could be attained with production helicopters. In the past, such specifications had to be 
modified considerably to a level realizable with current technology. In 1976, the NASA Research and 
Technology Advisory Council Subpanel on Helicopter Technology recommended a cabin vibration level 
of 0.02 g. This level could be achieved only if there were a major technological breakthrough in 
vibration control methods. 

Recently, a major advance in vibration control technology appears possible with the development 
of higher harmonic control (HHC), also known as multicyclic control. In contrast to passive vibration 
control systems which suppress the vibratory loads after they have been generated, an HHC system 
suppresses the vibrations at the source — the higher harmonic blade airloads. Higher harmonic control 
is an active vibration control system wherein the vibratory airloads are controlled at the higher har- 
monics of rotor rotational speed. Operated through an arrangement of sensors and microprocessors, 
higher harmonic inputs are applied to the blade generating new unsteady aerodynamic loads which, in 
combination with the resultant inertial loads, cancel the existing vibratory blade loads that cause the 
airframe vibration. 

There are many ways to implement HHC on a rotor system. The most popular approach is blade 
root feathering using swashplate oscillations. Using servo-actuators, the swashplate is excited in the 
collective, longitudinal cyclic, and lateral cyclic modes at per rev, resulting in blade pitch oscillations 
at three distinct frequencies of N^ — 1, N^, and A^ + l per rev in the blade rotating frame. For other HHC 
systems, the higher harmonic inputs are applied using jet flaps (refs. 4-7), partial-span flaps (refs. 8-11), 
circulation control rotor (refs. 12 and 13), and individual blade control system (refs. 14 and 15). 
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1.2 Literature Survey on HHC 


The concept of higher harmonic control was first proposed in the late forties, but it was not until 
1961 when it was conceived as a potential vibration control method for helicopters. In 1948, Winson 
(ref. 16) suggested that the higher harmonic blade pitch could be applied to rotors tested on whirl towers 
to simulate the oscillatory blade stresses induced in forward flight. This idea was pursued by Jensen 
(ref. 17) in 1951, who applied second harmonic control to the S-52 metal main rotor blades operated 
on the Sikorsky rotor test stand. 

In 1952, Stewart (ref. 18) suggested that the second harmonic control applied to rotors in forward 
flight could be used to delay the onset of retreating blade stall. Based on his analysis, Stewart discovered 
that the rotor disk loading could efficiently be redistributed using second harmonic control. For a 
particular flight condition, such loading redistribution could be adjusted to avoid retreating blade stall. 
The resulting effect would be to raise the speed limitation of helicopters. According to his analysis, 
helicopter speed limit could be increased by the order of 0.1 in terms of advance ratio. His analysis 
was, however, limited to a rigid flapping blade, and the airloads were calculated with quasi-steady 
aerodynamics and a uniform inflow distribution. Therefore, the severe aerodynamic effects associated 
with the transonic effects, separated flow conditions, and the participation of blade elastic modes in the 
rotor aeroelastic responses were all neglected. 

In 1961, Arcidiacono (ref. 19) extended Stewart’s analysis by considering both second and higher 
harmonic blade pitch control. In his report, Arcidiacono noted that a combination of second and third 
harmonic blade feathering could be used to delay the onset of retreating blade stall to an even higher 
advance ratio than that reported by Stewart. Like Stewart, Arcidiacono did not consider the effects of 
higher harmonic inputs on the rotor vibratory hub loads. 

In 1961, for the first time, a flight test program was carried out to investigate the feasibility of using 
HHC for vibration reduction on an UH-1A two-bladed teetering rotor helicopter (ref. 20). Prompted by 
some earlier unpublished analyses and a design of the simple rotor head mechanism to generate second 
harmonic blade pitch, Bell Helicopter Company conducted a series of flight tests to determine the effects 
of HHC on rotor performance, blade airloads, blade and control loads, hub loads, and airframe vibration. 

Investigating HHC effects on rotor performance, these investigators noted that no reduction in 
shaft torque was observed. Further investigation indicated that the drag reduction in the retreating side 
was accompanied by an increase in profile drag in the fore and aft portions of the rotor disk. These 
results confirmed Stewart’s prediction that the second higher harmonic control could be used to alter 
the rotor disk loading, although no reduction in power was observed at 100 knots. The reduction in the 
stall-induced power was later found to be commensurate with an increase in the profile power due to 
changes in the blade section angle of attack. 

On vibration, these authors stated that none of the combinations of second harmonic control am- 
plitude and phase could completely suppress the vibration at any one of three fuselage stations being 
monitored. At best, the vibration measured at the station located right beneath the rotor hub was re- 
duced by about 50 percent. These authors pointed out that the complex interactions between the various 
components of the vibratory hub loads were the cause of incomplete vibration reduction and concluded 
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that more comprehensive understanding of HHC effects was required before the system could become 
fully effective. 

Since then, several investigators had conducted analytical simulation of HHC for vibration reduc- 
tion. Daughaday (ref. 21) presented the analysis using a type of higher harmonic blade pitch control 
which required both temporal and spanwise variations. These higher harmonic blade pitches were 
necessary to reduce a set of even harmonic shaftwise hub forces for a two-bladed teetering rotor. In 
particular, this type of control schedule is quite impractical because it requires the adjustability of the 
blade pitch angle along the blade span. In this analysis, the rotor induced inflow is calculated using a 
vortex wake analysis, and the effects of shed vorticity were also included. In an analysis similar to that 
of Daughaday, Shaw (ref. 22) carried out an HHC investigation including the effects of third harmonic 
(rigid) blade root feathering on the vibratory hub forces for a three-bladed rotor. 

Both of these results indicated that the higher harmonic shaftwise root forces could be eliminated 
using HHC and that the load suppressions were caused by the mutual cancellation of the higher har- 
monic aerodynamic and inertial components. Furthermore, interharmonic couplings — the interaction 
between the harmonics of inputs with another harmonics of outputs — were observed by both authors. 
Interharmonic couplings are expected to exist for all rotor systems. For example, even with the simplest 
mathematical model for a rigid flapping blade rotor in forward flight, the blade motion is described by a 
linear ordinary differential equation with periodically time-varying coefficients. For such an equation, it 
is well established that the input (or forcing) at one harmonic generates a multiharmonic output. Finally, 
Daughaday noted in his report that the HHC inputs required for vibratory hub load suppression changed 
significantly with forward speed. 

Daughaday’s analysis on the two-bladed teetering rotor was later extended by Balcerak and Erikson 
(ref. 23). The modifications included the torsional degree of freedom, which had been previously 
neglected by Daughaday, and the calculation of the blade inplane root shears. These authors calculated 
the higher harmonic blade pitch required to eliminate a set of even harmonics of shaftwise hub force 
together with a set of odd harmonics of inplane hub shear. Results from this study pointed out the 
importance of including the torsional degree of freedom, especially when the blade fundamental torsion 
frequency was close to the blade pitch harmonics. Also, the rotor aerodynamic performance was found 
to be relatively unaffected by the HHC inputs. 

In 1972, the first wind tunnel investigation of HHC was reported by McCloud (refs. 4 and 5). 
The model rotor was a two-bladed teetering rotor with propulsive jet flaps. This rotor used a large jet 
flow expelled from the blade trailing edge to propel the rotor and applied the HHC through the angular 
deflection of the jet flow. In the model rotor, the jet flap was extended over the outer 30 percent of the 
blade radius. In the experiment, individual higher harmonic input components were applied to determine 
a set of influence coefficients relating the vibratory blade and hub load harmonics to the input harmonics. 
This set of influence coefficients formed a transfer matrix, which was then used to determine the higher 
harmonic input magnitudes and phases required to minimize the sum of the squares of the vibratory 
loads. These results showed that, in general, the vibratory hub load reduction was accompanied by 
an increase in the blade bending moments and vice versa. As noted by earlier investigators, the HHC 
inputs required for the suppression of vibratory hub loads were found to vary appreciably with the rotor 
forward speed. These experimental results were later supported by the analytical studies performed by 
Piziali and Trenka (ref. 6) and Trenka (ref. 7). 
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Sissingh and Donham also performed an HHC investigation in a wind tunnel using a four-bladed 
hingeless rotor (ref. 24). Due to the actuator hardware constraints, the rotor tip speed was significantly 
low compared to a full-scale helicopter. The resulting blade vibration frequencies were therefore very 
different from those of full-scale hingeless rotors. In the wind tunnel experiment, the higher harmonic 
blade root feathering was applied using swashplate excitations. The nonrotating swashplate was excited 
with servo-actuators oscillated at the fourth harmonic of the rotor rotational speed. This results in 
blade pitch harmonics at 3, 4, and 5 per rev which were superimposed on the basic collective and 
cyclic controls to suppress the third and fifth harmonics of the blade flap bending moments measured 
near the blade root. In effect, this suppression scheme tried to null out the fourth harmonics of the 
hub pitching and rolling moments in the fixed system. The vibratory hub loads, however, were not 
measured. The HHC model used for vibratory hub load suppression was similar to the one used by 
McCloud (ref. 5). With the applied HHC inputs, the results showed that these vibratory flap bending 
moments were drastically reduced and the HHC magnitudes generally increased with increasing advance 
ratios. 

Starting in 1975, a series of wind tunnel investigations of HHC was conducted by McHugh and 
Shaw (refs. 25 and 26). The four-bladed soft-inplane hingeless rotor model being tested was fully instru- 
mented so that secondary effects of HHC could be investigated. HHC inputs were applied by swashplate 
oscillations with servo-actuators. The vibratory hub loads were measured using a strain-gauge balance 
located beneath the actuators, in the fixed frame. Unlike Sissingh and Donham’s experiment, the later 
experiments of McHugh and Shaw were conducted at close to the full-scale tip speed. Results from 
the experiment showed that the vibratory hub moments could be suppressed effectively with no signifi- 
cant increase in blade stresses. Also, all five components of the fourth harmonics hub loads — vertical, 
lateral, and longitudinal forces and pitching and rolling moments — could be reduced simultaneously 
with the three HHC inputs. However, the vibratory hub shears could not be completely suppressed for 
certain flight conditions, and this was attributed to the amplitude limitations of the servo- actuators. No 
significant penalties in rotor performance were observed in the HHC experiment. Noting that the HHC 
amplitudes and phases varied appreciably with forward speed, the authors suggested that an adaptive 
mechanism would be necessary to adjust the HHC inputs as required for varying flight conditions. 

In 1975, McCloud presented an analytical investigation of HHC on a four-bladed articulated 
controllable-twist rotor (ref. 8). The analysis used an aeroelastic program developed for controllable- 
twist rotors by Kaman Aerospace. McCloud used the analysis to obtain the influence coefficients, 
relating the vibratory blade and hub load harmonics to the higher harmonic inputs, for the HHC model. 
In the study, McCloud confirmed the validity of the linear assumption made on the HHC input-output 
model. Using a weighting approach, he further showed that the suppression of the vibratory hub load 
could be accompanied by a 50 percent reduction in the blade bending stresses using the first four har- 
monics of the servo-flap inputs. Another analytical investigation of HHC on the controllable-twist rotor 
was presented by Lemnios and Dunn (ref. 9). 

McCloud’s HHC investigation on the four-bladed articulated controllable-twist rotor was carried 
out with a full-scale test in the Ames 40- by 80-Foot Wind Tunnel (ref. 10). In this experiment, the 
rotor blades were fully instrumented, and the rotor tip speed was about one-fourth that of production 
helicopters. The vibratory blade loads were obtained for different higher harmonic servo-flap deflections 
and at varying flight conditions, and the data were stored for a post-test analysis. In the data analysis, 
model reductions were performed using a regression analysis. The HHC model was then implemented 
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in a simulation program to compute the vibratory blade load reductions. The HHC mathematical model 
and controllers used in the simulation were identical to those of McCloud’s analytical study presented 
in reference 8. Based on the simulated results, the authors reported that significant reductions occurred 
in both the blade flapwise bending moment and the control system loads using HHC. In addition, the 
data analysis also indicated that rotor performance could be improved using HHC. 

In another HHC investigation of the same rotor, reported by McCloud and Brown in 1980 (ref. 11), 
the simulation included fuselage vibration measurements. In the earlier test, vibrations were measured 
using accelerators mounted on the rotor test apparatus (RTA). While varying the weighting options, they 
discovered that significant vibration reduction could be achieved with HHC at selected locations on the 
RTA. These authors also noted that components of the transfer matrix varied appreciably with forward 
speed and, like Shaw, suggested that a real-time parameter identification scheme would be required for 
practical HHC application. 

In a paper reported by Hammond (ref. 27) and another one presented by Wood, Powers, and 
Hammond (ref. 28), results from the wind tunnel investigation of HHC on a four-bladed articulated 
model rotor were described. Like previous wind tunnel investigations, the HHC inputs were applied 
using swashplate oscillations. The vibratory hub loads were measured by a six-component balance 
mounted below the rotor hub. Using several trial inputs, these investigators were able to compute the 
fourth harmonic blade pitch required to suppress completely the 4 per rev vertical hub force. Then the 
authors outlined a general algorithm to acquire and control the multi-input multi-output HHC model 
and discussed the practical implementation of HHC in flight tests. In fact, the practical design issues 
concerned with the application of HHC to a production helicopter were later investigated by Powers 
(ref. 29). 

The potential of HHC for future helicopter developments was reviewed by McCloud in 1980 
(ref. 30). The potential benefits of HHC as an active control system for helicopter vibration suppres- 
sion, blade stress reduction, and rotor performance improvement were outlined. Several types of rotors 
were reviewed as potential candidates for HHC applications including conventional, controllable-twist, 
propulsive and nonpropulsive jet flaps, circulation control, and advanced-controllable-twist rotors. Fur- 
thermore, practical implementations of HHC in flight test were also discussed. 

In 1980, Shaw presented the results of a comprehensive analytical investigation of HHC comparing 
the potential benefits of servo-flap versus conventional blade root feathering, and studied the development 
of an automatic in-flight adaptive algorithm (ref. 31). The investigation was based on a coupled modal 
analysis and included a vortex wake induced inflow calculation. The shed wake effects were included 
using an approximation to the Theodorsen lift deficiency function. In the modal analysis, the inplane 
degree of freedom was uncoupled from the out-of-plane and the torsional degrees of freedom. The 
transfer matrix approach was adopted to relate the higher harmonic hub loads to the HHC inputs. 

The analytical results indicated that nonlinearities in the HHC input-output model were small, on 
the order of 10 percent. Also, Shaw showed that significant reduction in vibration could be achieved 
with small to moderate increases in blade stresses. Investigating the distributions of the baseline and 
the incremental vibratory blade loads induced by HHC, Shaw pointed out that the vibration suppression 
was caused by the mutual cancellation between the aerodynamic and the inertial components of the 
transmitted vibratory loads at the blade roots. Cancellation of the blade root loads by using blade root 


8 



feathering generally caused an increase in both the aerodynamic and inertial loading distribution from 
the baseline values, thus resulting in an increase in blade stresses. Further investigation indicated that 
the increases in blade stress could be avoided if partial span flap, as opposed to blade root feathering, 
was used for vibration suppression. 

Due to the HHC inputs, it was also found by Shaw that the control loads were increased by roughly 
30 percent, and the change in rotor performance was negligible. For changing flight conditions, the 
closed-loop controllers with fixed gain performed satisfactorily over an advance ratio range as wide as 
0.2. For the cases when this controller performed poorly, an adaptive gain controller was used. For the 
adaptive controller, the model parameters were estimated using the Kalman filter. Simulation results 
showed that the adaptive controller performed very well for varying flight conditions. Furthermore, the 
predicted actuator power requirement for HHC was less than 1 percent of the rotor induced power. 

Starting in 1979, Shaw and Albion conducted wind tunnel investigation of HHC on a four-bladed 
hingeless rotor model in the Boeing V/STOL Wind Tunnel (refs. 32 and 33). The 10 foot diameter rotor 
was Mach scaled and operated at the full-scale tip speed. The HHC inputs were applied through the 
fixed-system swashplate excitation using four equally spaced hydraulic servo-actuators. The closed-loop 
controller suppressed simultaneously up to 90 percent the 4 per rev vertical hub shear, hub pitching, 
and rolling moments. The effectiveness of HHC was demonstrated for the advance ratios ranging from 
0.1 to 0.3. The closed-loop system transient behaviors were investigated by applying a step disturbance 
through a swashplate command. Results indicated that the disturbance was suppressed within two rotor 
revolutions, or approximately 0.1 second. These results confirmed the quasi-static assumption made on 
the HHC model. Also, it was noted that the rotor trim states were not affected by the HHC inputs. 
However, the maximum penalty in rotor performance was large, on the order of 4 percent. 

In 1980, Yen presented an analytical investigation of HHC on a two-bladed teetering and a four- 
bladed hingeless rotor (ref. 34). The analysis was based on C81, the Global Flight Simulation Program. 
For the two-bladed teetering rotor, trial amplitudes and phases of the second harmonic inputs were 
applied to suppress the second harmonics vertical hub shear. At the advance ratio of 0.3, the main rotor 
shaft torque was increased due to the applied HHC by 1 1 percent. The increase in shaft torque was 
attributed to blade stall which occurred in the fourth quadrant of the rotor disk. 

For the four-bladed hingeless rotor, third, fourth, and fifth harmonics of blade pitch were simulta- 
neously applied to suppress the vibratory hub loads. This resulted in reduction of more than 90 percent 
in the 4 per rev vertical hub shear, hub pitching, and hub rolling moments. For this hingeless rotor, 
the combined HHC inputs required were about 2.5 degrees, and the HHC penalty in rotor performance 
was found to be small. These results indicated that the second harmonic blade pitch, implemented on 
a two-bladed rotor, changed the rotor loading distribution and resulted in blade stall. Such effects had 
previously been suggested by Stewart (ref. 18). Furthermore, these results also implied that the blade 
pitch harmonics above the second, used for the vibratory hub load suppression, had a small or negligible 
effect on rotor performance. 

In 1980, Hammond (ref. 35) reported a wind tunnel investigation of HHC on a dynamically Mach- 
scaled four-bladed articulated model rotor. The model rotor tip Mach number matched that of a full-scale 
value. As in Hammond’s previous wind tunnel experiments, the HHC inputs were applied using swash- 
plate oscillations. The transfer matrix approach was used for HHC modeling. Self-tuning regulators 
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were used for vibration suppression, and the model parameters were estimated using the Kalman fil- 
tering techniques. For a range of advance ratios between 0.2 to 0.4, the 4 per rev vertical hub shear 
and hub pitching moment were suppressed up to 70 percent from the baseline values. In general, the 
4 per rev hub rolling moment could not be suppressed, however. Like other investigators, Hammond 
also reported a moderate increase in blade loads due to HHC. 

Hammond and Hollenbaugh (ref. 36) investigated the application of different HHC control schemes 
to the vibration suppression problem. The HHC control schemes included the application of a single- 
input, multi-input, and multi-input with on-line adjustable blade twist. The last control scheme was quite 
similar to the one proposed by Daughaday (ref. 21) in which the added blade pitch required both spatial 
and temporal variations. The vibration suppression problem was solved using a nonlinear programming 
algorithm in which the quality criteria to be minimized included the sum of the peak-to-peak values of 
hub loads. The mathematical model used for the vibration problem was quite simple. It used either 
a rigid flapping or a coupled elastic flap-lag blade model, quasi-steady aerodynamics, and a linear 
inflow model. The results showed that for suppression of the oscillatory hub forces and moments on 
a four-bladed hingeless rotor, simultaneous application of the third, fourth, and fifth harmonics of the 
rigid blade pitch yielded good results. Better results were obtained with the addition of the variable 
time-dependent blade twist. However, a complete suppression of the oscillatory hub loads could not be 
achieved. 

In 1981, Molusis, Hammond, and Cline (ref. 37) examined several HHC algorithms for vibration 
suppression. The performance of these algorithms was evaluated through wind tunnel testing. The 
model rotor and experimental setup were identical to that of Hammond’s previous wind tunnel test 
(ref. 35). The HHC controllers were configured to suppress the fourth harmonic vertical, longitudinal, 
and lateral signals from a triaxial accelerator mounted beneath the rotor in the fixed frame. The rotor 
advance ratio was varied between 0.2 to 0.4, and the HHC system was modeled using the transfer matrix 
approach. 

The HHC controllers tested could be separated into two types: the adaptive controllers and the gain 
schedule controllers. For the adaptive controllers, both the baseline vibration and the transfer matrix 
were estimated using Kalman filtering techniques. These adaptive controllers could be classified further 
depending on how the model parameters were updated. The adaptive controllers being evaluated were 
the deterministic and the cautious controllers. In particular, the deterministic (adaptive) controller used 
the estimated model parameters directly. For the cautious controller, the model parameters were also 
estimated using the Kalman filter. However, the cautious terms, which were calculated based on the 
covariances of the estimated parameters, were included in the control algorithm. If the parameter estima- 
tion was uncertain, then the controller decreased the gain in the control function. For the scheduled-gain 
controller, the transfer matrices were assumed to be known and were calculated off line. The two types 
of gain-scheduled controller being tested for different flight conditions in the wind tunnel included the 
perturbation and the proportional controllers. For the perturbation controller, the HHC model was linear 
about a nominal control value, whereas, for the proportional controller, the transfer matrix was fixed for 
the complete range of control values. 

Wind tunnel results indicated that both gain-scheduled controllers performed poorly, and there was 
no definite reduction in the vibration level. It was suggested that a nonlinear behavior might be the 
cause for such poor performance. For the deterministic (adaptive) controller, the steady state vibration 
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level was reduced significantly. However, large response overshoots occurred before the vibration 
converged to the steady state value. The cautious controller offered the best performance among the 
four controllers. After initialization, the vibration converged smoothly to a minimum value. Over the 
range of forward speed tested, the cautious controller was able to suppress up to 90 percent of the fourth 
harmonic vertical and longitudinal accelerations. However, the 4 per rev lateral vibratory component 
could not be suppressed, in general. 

In an analytical study, Molusis (refs. 38 and 39) investigated the inability of the HHC controllers to 
suppress the 4 per rev hub rolling moment as reported in reference 35 and the 4 per rev lateral acceleration 
as discussed in reference 37. Using the G400 simulation program for aeroelastic analysis (ref. 40), 
Molusis claimed that nonlinearity in the HHC modeling was a major factor for the incomplete vibration 
suppression. The nonlinearity was found to be more severe at the high forward speed. In fact, the 
incomplete vibration suppression was attributed to two factors: (1) the convergences of HHC algorithms 
to local minima, and (2) the constraints on the HHC amplitudes (the combined HHC amplitudes in 
Hammond s experiment were limited to 1.5 degrees). Furthermore, Molusis also noted that the parameter 
identification scheme using Kalman filtering technique could, in certain cases, lead to divergent HHC 
transient responses. 


In 1982, Wood et al. (refs. 41-43) conducted an HHC flight test program on a modified OH-6A 
helicopter. The OH-6A is a light helicopter with a gross weight of about 3,000 lb. The flight test 
program was supported by the earlier wind tunnel results performed by Hammond. The HHC inputs 
were implemented by blade root feathering using swashplate oscillations, and a triaxial accelerator 
was mounted beneath the pilot seat to sense and feed back the fourth harmonic vertical, longitudinal, 
and lateral vibrations to the controller. The flight speed was varied between zero to 100 knots, while 
HHC was applied using both open-loop (manual) and closed-loop controllers. For the closed-loop 
operation, the cautious controller presented in reference 37 was employed. Flight test results confirmed 
the effectiveness of HHC, as had been observed in wind tunnels. Up to a 90 percent reduction in 
vibration was obtained with combined HHC amplitudes less than 1 degree. The blade and control 
system fatigue loads were increased due to HHC; however, the load penalties were within the design 
limits. Shaft torque measurement showed no adverse effects of HHC on the aircraft performance. 

Johnson (ref. 44) reviewed the class of HHC algorithms for vibration and blade load reduction 
proposed to date. As in previous investigations, a linear, quasi-static, frequency domain transfer matrix 
was adopted for the HHC modeling. Furthermore, the HHC models were classified into global (linear) 
and local (piecewise linear) models. Both open-loop and closed-loop controllers were considered. 
The criteria for vibration suppression or blade load reduction were expressed as a minimum variance 
or a quadratic performance index. The controllers were further classified according to the parameter 
identification scheme being used. As in previous reports, these included the deterministic, gain-updated, 
cautious, and dual controllers. 

The performance of the controllers presented by Johnson was evaluated analytically by Chopra 
and McCloud (ref. 45). The performance evaluation was based on a numerical simulation which used 
the data obtained from McCloud’s earlier wind tunnel results (ref. 11). They studied the effects of 
measurement noise, initial estimated error in the model parameters, control rate constraint, and variation 
in flight speed on the performance of four different controllers. In particular, these controllers were 
open loop and closed loop with off-line and on-line identification. The on-line identification scheme 
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was performed using the Kalman filter. Both global and local models were simulated. The closed-loop 
controller with a local model using on-line identification technique had the best performance. However, 
for the cases with large initial errors in the transfer matrix, large overshoots were found in the transient 
responses using this controller. For such cases, these authors proposed two different startup procedures 
to improve the transient responses: (1) trial HHC input sequences would be injected into the system 
to compute a better estimate of the transfer matrix using either the finite-difference or the least-square 
methods, and (2) the control rate constraint value would be made initially large to avoid the undesirable 
control excursion and eventually reduced for improved convergence characteristics. 

In 1985, a successful wind tunnel demonstration of HHC was performed by Shaw et al. (ref. 46). 
The experimental investigation was conducted using a three-bladed dynamically scaled model rotor of 
the CH-47D helicopter tested in the Boeing Helicopters Wind Tunnel. The model rotor operated at 
full-scale tip speed, and the vibratory hub loads were measured using a strain-gauge balance mounted 
in the rotating frame. The higher harmonic blade feathering at the second, third, and fourth harmonics 
was applied to suppress the 3 per rev vertical hub force and the 2 and 4 per rev rotating inplane hub 
shears. Open-loop tests were conducted to obtain the transfer matrix values at several flight conditions. 
These transfer matrices were stored and later used with the fixed- or scheduled-gain controllers. 

Wind tunnel results showed that HHC produced more than 90 percent reduction of all three vibratory 
hub shear components. Such reductions were obtained using the fixed-gain controller with a local model, 
in which the transfer matrices were calculated off line. Vibration reduction was shown for the complete 
range of flight conditions being tested, up to a forward speed of 188 knots. The suppression of the 
vibratory hub loads was accompanied by small to moderate changes in the blade and control system 
fatigue loads as well as a negligible change in rotor performance. However, the combined HHC 
amplitudes required for vibratory load reduction were relatively large — on the order of 2 and 3 degrees, 
respectively, at the forward speeds of 135 and 188 knots. These HHC amplitudes were significantly 
larger than those reported in reference 41. The reason, as suggested by Shaw, was due to the differences 
in blade dynamics between the two rotors. 

Testing different HHC algorithms, these authors found that the scheduled-gain controller performed 
as well as the fixed-gain controller at varying flight conditions. It was interesting since the scheduled- 
gain controller selected the transfer matrices automatically from a data base, while the transfer matrices 
were selected manually for the fixed-gain controller. Also, both controllers used the local model. The 
adaptive controllers, similar to those reported in reference 37, were also tested; nevertheless, they 
generally exhibited poor performance. For certain fixed flight conditions, the adaptive controller was 
successful in suppressing the vibratory hub loads. For changing flight conditions, however, the transfer 
matrices identified by this controller diverged for certain values of the covariances of the measurement 
and process noises. Furthermore, some initial testing of the global adaptive controller indicated that it 
was either unstable or ineffective in suppressing the vibratory loads. 

The HHC effectiveness shown in this wind tunnel experiment had other implications. First, it 
demonstrated that the vibratory hub loads could be suppressed with a relatively simple HHC controller — 
the fixed-gain controller. Second, the nonlinearities in the HHC model, if these existed, did not prevent 
the near complete suppression of the vibratory hub loads. In fact, open-loop results indicated that the 
vibratory responses were linear with HHC amplitudes, up to 3 degrees. These results were at variance 
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with those reported in references 37 and 38, and the differences may be attributed to the actual hardware 
implementation, the servo-actuator accuracy, and the data processing techniques. 

In 1985, Lehmann reported an HHC investigation conducted in the the Deutsch-Niederlandischer 
Windkanal (DNW) (ref. 47). This wind tunnel investigation was supported by an earlier analytical study 
performed by Jacob and Lehmann (ref. 48). The higher harmonic blade pitch was applied to suppress 
the vibratory hub loads on a four-bladed hingeless model rotor of the BO- 105 helicopter. The load 
balance was mounted in the fixed frame, beneath the rotor. The airspeed was varied between 10 and 
50 meters per second which corresponded to advance ratios of 0.04 to 0.22. The HHC results indicated 
that all five vibratory hub load components— 4 per rev vertical, longitudinal, and lateral hub forces and 
4 per rev pitching and rolling moments — could be reduced simultaneously but could not be completely 
suppressed. Furthermore, Lehmann found that the rotor trim states were affected by the applied HHC 
inputs. 

Analytical investigation of HHC for suppression of fuselage vibration was also carried out by Taylor, 
Farrar, and Miao (ref. 49). The HHC system being considered also employed swashplate oscillations. 
A four-bladed single rotor helicopter was simulated using the G400 rotor aeroelastic analysis. The 
controllers used were similar to those described in reference 35. Simulated results indicated that HHC 
produced up to 80 to 90 percent reduction in fuselage vibration. Secondary effects of HHC included a 
small reduction in rotor performance— on the order of 1 to 3 percent increase in rotor shaft torque— and 
small to moderate increase in blade stresses. The influences of measurement noise included in the 
simulation indicated that a 15 percent noise-to-signal ratio was acceptable for the HHC design. With 
the rotor system considered, the predicted HHC inputs required were on the order of 3 degrees. 

Another comprehensive HHC simulation, also using the G400 program, was carried out by Davis 
(ref. 50). This study aimed at comparing the performance of three different controllers: deterministic, 
cautious, and dual. Both local and global models were used. For steady flight conditions, Davis noted 
that, with proper tuning, there was no distinguishable difference between the performance of these three 
controllers. The vibration reductions were 75 to 95 percent from the baseline values, and these were 
achieved with less than 1 degree of the combined HHC inputs. For transient maneuvers, it was found 
that the global controllers generally performed better than the controllers using the local model. 

The effects of nonlinearity and interharmonic coupling of the HHC system were also investigated 
in reference 50. The major source of nonlinearity in the HHC model came from the blade bending 
mode responses, while the interharmonic coupling arose from both the blade responses and the higher 
harmonic airloads. Both of these effects might degrade the controllers’ performance in terms of the 
HHC amplitude required and the vibration convergent rate. By selecting proper weighting factors on the 
HHC inputs at different flight conditions, the author noted that such degrading effects might be reduced. 

The above simulation led to a study for implementing HHC on the Sikorsky S-76 helicopter 
(ref. 51). Analytical studies, design analysis, and risk reduction tests were all performed prior to the 
full-scale flight demonstration. The flight test program was designed to demonstrate HHC effectiveness 
on a modem medium-weight helicopter. 
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In 1986, flight demonstration of HHC on the S-76 helicopter was conducted and described in 
references 52 and 53. For flight tests, passive vibration control devices were either removed or switched 
off. Swashplate excitation was used to activate the higher harmonic blade pitch, and the cabin vibration 
level was monitored using accelerators mounted at several locations in the aircraft. Using open-loop 
testing, this was the first flight test to demonstrate the potential of HHC on a 10,000 lb helicopter for 
a range of forward speed up to 150 knots. With HHC, significant reductions in the fuselage vibration 
were attained. The effects on aircraft performance were found to be negligible, while the penalties 
in blade and control system fatigue loads were moderate and remained well within the design limits. 
Due to actuator power limitations, however, only a small percentage in cabin vibration reduction was 
achieved at forward speeds above 100 knots. 

Also in 1986, Polychroniadis and Achache (refs. 54 and 55) presented an HHC flight test on an 
SA 349 Gazelle experimental research aircraft, a lightweight (4,500 lb) helicopter. The HHC inputs 
were implemented using swashplate excitation. Self-adaptive controllers used feedback signals from the 
vibration sensors placed at selected locations in the aircraft cabin. The controllers used in these flight 
tests were similar to the adaptive controllers presented in reference 37. For the flight tests, passive 
vibration control devices were removed from the aircraft. With HHC, a 70 to 90 percent reduction in 
vibration was achieved at forward speeds up to 250 km per hour (or 135 knots). These results indicated 
a superior performance in terms of vibration suppression using HHC over the previously used passive 
vibration control devices. The interactions between the vibratory hub loads and the airframe vibration 
were also investigated. These authors found that the vibration suppression using HHC was achieved by 
significant reductions in the second harmonic rotating inplane hub forces and moments. At the same 
time, the fourth harmonic components of the inplane hub loads were increased slightly. 

In 1988, Jacklin (ref. 56) presented a numerical simulation evaluating the performance of five 
system identification techniques for HHC. These included the weighted least square method, the classical 
Kalman filter, a generalized Kalman filter, the classical least mean square filter, and the generalized least 
mean square filter. Both global and local models were considered. Results from both open-loop and 
closed-loop simulations were shown. The algorithms were evaluated in terms of accuracy, stability, 
convergence properties, computational speed, relative difficulties for practical implementation, and the 
level of vibration reduction. Noises and modeling errors were included in the simulation study. Based 
on the results, Jacklin concluded that the generalized least mean square method using the local model 
in the closed-loop operation was the most promising identification method for HHC. 

In 1988, Sangha (ref. 57) presented an analytical investigation of HHC using the Rotor/Airframe 
Comprehensive Aeroelastic Program (RACAP). This study concentrated on the influences of the open- 
loop HHC on the fuselage vibration. Based on a correlation study, Sangha pointed out the importance of 
including hub motions for the prediction of fuselage vibration. In the design of an HHC system, he noted 
that it was the modal properties of the fuselage that were important, not those of the blade. However, 
such conclusion was made based only on the parametric study of the control system stiffness on the 
4 per rev vibration at the pilot seat using the open-loop controller. However, other blade characteristics 
such as distributions of blade inertia and stiffness, which have strong influences on the HHC responses, 
were not addressed. Also, Sangha showed that the vibration reduction was caused primarily by the 
phase variations rather than the magnitude variations of the vibratory hub loads. 
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In 1989, Robinson and Friedmann (ref. 58) conducted an analytical investigation of HHC on 
a generic four-bladed hingeless rotor. The results were calculated using the coupled flap-lag-torsion 
aeroelastic analysis developed at the University of California, Los Angeles. Both quasi-steady and 
unsteady aerodynamic results were presented (ref. 59). The unsteady aerodynamic model employed was 
a finite-state time-domain version of the approximated Wagner or Theodorsen lift deficiency functions 
(refs. 60 and 61). These authors showed significant differences between the quasi-steady aerodynamic 
and the unsteady aerodynamic results for the computed 4 per rev vertical hub shear, and the differences 
were more pronounced at the high advance ratios. The results for the inplane hub shears were essentially 
identical between these two aerodynamic models. 

Deterministic and cautious controllers were used to suppress the predicted vibratory hub loads 
using both the global and the local models. From the simulated results, no significant differences were 
found among the performance of the controllers being tested. The computed HHC inputs for vibration 
suppression were found to be noticeably different depending on which aerodynamic model was being 
used. 

Also in 1989, Welsh and Blackwell presented an HHC study on the X-Wing (four-bladed) model 
rotor conducted in the Large Subsonic Wind Tunnel at the United Technology Research Center (ref. 1 2). 
For the X-Wing rotor, both trim controls and higher harmonic inputs were applied using the circulation 
control concept. In a circulation control airfoil, the combination of trailing edge blowing and Coanda 
effect can create more lift than conventional airfoils. In the model rotor, compressed air was supplied 
to a stationary annular plenum and passed through a series of stationary valves to receivers feeding 
the leading edge and trailing edge ducts mounted on the blade. The ducts in the blade carried the 
compressed air down the span and vented the air through slots to the Coanda surface of the airfoil. 
Unlike the standard swashplate excitation with mechanical linkages and actuators, the higher harmonic 
inputs could be applied without any moving part and at any harmonic of the rotor rotational speed. For 
this investigation, HHC was applied using trailing edge blowing at 2, 3, 4, and 5 per rev. The vibration 
was measured using a three-axis accelerometer mounted on the stand. Also, a balance was used to 
measure the vibratory hub loads, and the blade were fully instrumented with strain gauges so that the 
blade stress could be monitored. Reductions in both the vibration and the blade stresses were achieved 
using a fixed-gain controller in which the transfer matrix was calculated off line. The forward speed 
was varied between 60 to 175 knots, and the rotor operated at 50 and 100 percent full-scale tip speed. 

With the HHC, the vibration was suppressed between 70 and 95 percent from the baseline level, and 
the blade loads were also reduced simultaneously. The blade load reduction was primarily caused by the 
second harmonic control, typically not realizable on a four-bladed rotor using conventional swashplate 
oscillation. Also, it was found that the coupling between the trim control and the HHC was so severe 
that such effects had to be included into the control algorithm to maintain the rotor in trim. Analytical 
investigation of HHC on the X-Wing circulation control rotor was also performed by Abramson (ref. 1 3). 

Up to this point, all the active vibration control algorithms discussed are in the frequency domain. 
In 1981, Du Val, Gregory, and Gupta (ref. 62) advocated a time domain approach for the active control 
of helicopter vibration. In fact, this method was an extension of the linear quadratic Gaussian (LQG) 
regulator with a frequency shaped cost function. In essence, the proposed approach used a constant 
regulator with a filter in the feedback loop. The feedback variables, such as the fuselage accelerations, 
were passed through an undamped second-order oscillator, tuned to the frequency of TV/, per rev. At 
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the filter resonant frequency, the regulator locked onto the phase and magnitude of the feedback signal 
without the need for harmonic analysis. The proposed technique also offered other advantages over 
the frequency domain approach. Since state feedback was used, on-line identification of the model 
parameters would not be necessary. Furthermore, linear system theory techniques could be employed 
for the design of a stable system. 

Using a simple simulation program for the vibration control model, these investigators were able to 
show that this control scheme performed extremely well. With the feedback controllers generating the 
collective, longitudinal, and lateral cyclics and tail rotor controls, the reduction in vibration levels was 
more than 90 percent. However, the controller used for vibration suppression was initially designed for 
the hovering condition, and, as such, the effects of the periodic-time-varying coefficients in the plant 
model had been neglected. This control scheme was also used to suppress the vibration on an X-Wing 
model rotor (ref. 63). 

In 1989, Wereley and Hall (ref. 64) presented an interesting and alternative view of the active 
vibration control algorithms currently being proposed. This study focused on evaluating the performance 
of these control algorithms within the context of classical control theory. These authors pointed out that 
all these active vibration algorithms, either in the time or the frequency domain, were in fact quite similar 
and could be considered, in the context of linear control theory, as the classical disturbance rejection 
problems. Furthermore, they claimed that the effectiveness of the HHC algorithms was rather limited 
by the quasi-static and the time-invariant assumptions imposed on the plant model. They also noted 
that all the HHC algorithms were fairly robust. And finally, they proposed two possible approaches for 
improving HHC performance. First, the periodic behavior of the rotor system should be included in the 
plant model. The second approach involved the development of a control algorithm for the particular 
periodic-time-varying plant model. However, these two approaches were still under investigation, and 
therefore no result was shown. 

Besides the works discussed above, significant contributions to HHC had been performed by Kretz 
(refs. 65—67). Application of a nonlinear programming method to HHC algorithm was discussed by 
Perley et al. (ref. 68). A feasibility study for the application of HHC to the advanced blade concept 
rotor was reported in reference reference 69. Prediction of HHC responses with a coupled rotor-airframe 
vibration control model was presented in reference 70. And finally, general discussions about HHC also 
appeared in references 1, 2, and 71-74. 


Besides the HHC using swashplate oscillation, several other vibration control methods using active 
control techniques were also investigated. These include individual blade control (IBC) (refs. 14 and 15) 
and active isolation (refs. 75 and 76) systems. In an IBC system, each individual blade is controlled by 
a separate servo-actuator, and, in fact, these servo-actuators are mounted on the rotating swashplate and 
thus rotate with the blades. The blade excitation can be achieved for any frequency. In reference 15, 
the IBC algorithms were formulated in the time domain using an extension of the standard linear 
quadratic regulator approach with the model-following systems. Without the restriction on the blade 
input frequency, the IBC system could be used to alleviate gust response and retreating blade stall, as 
well as to stabilize the rotor system. 

Active isolation is another attractive concept for suppressing helicopter vibration. In this system, 
servo-actuators are mounted in the airframe structures and configured to oscillate in such a way as 
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to cancel out the vibration. Either time or frequency domain control algorithms can be used for this 
system. If the frequency domain technique is used, then the control algorithm for this system is similar 
to that of the HHC. However, determination of the optimal locations and the number of servo-actuators 
required for actual implementation of active isolation is a central issue with this technique. Preliminary 
investigation using an analytical method has indicated that this vibration control system can perform 
virtually as well as the HHC system (ref. 71). 

1.3 Current Investigation 


Based on the discussion presented above, HHC can be considered to be one of the most effective 
vibration control techniques currently available for helicopters. Virtually all the investigations carried 
out to date on a number of different rotor systems have shown that HHC produced substantial reduction 
in either the vibratory hub loads or the fuselage vibration. Even though some of the HHC algorithms 
have been shown to perform better than others, each appears feasible for practical implementation. 

One problem which remains to be solved is the HHC amplitudes required for vibration suppression. 
While the results in reference 41 indicated that less than 1 degree in the combined HHC amplitudes was 
required, others — in particular, references 46 and 49 — have shown that approximately 3 degrees was 
needed to control the vibration on their rotors. This difference might be attributed to the differences in 
blade dynamics of these rotors. This problem has a practical significance since the HHC amplitudes are 
directly related to the power requirement of hydraulic servo-actuators. Indeed, for certain rotor designs, 
the weight penalty imposed by the actuation system can prevent the practical implementation of HHC 
on such rotors. On the other hand, an HHC system with a low amplitude requirement needs actuators 
of smaller weight and thus has a low weight penalty. Also, the maintenance cost is low since a small 
oscillation amplitude translates to a longer actuator seal life. Furthermore, it is desirable to minimize the 
secondary effects of the HHC. A large higher harmonic blade pitch may alter the helicopter trim state 
and thus affects the aircraft performance and handling qualities. Also, for aircraft operating near the 
flight boundary, HHC inputs with sufficiently large magnitude may induce a premature onset of blade 
stall. Since a successful HHC system must operate efficiently in a severe aerodynamic environment, its 
performance should be investigated over the whole flight envelope. 

The objectives of this study are to: 

1. Develop a comprehensive aeroelastic analysis to calculate the rotor responses in forward flight 
including consistent and complete unsteady aerodynamic and wake effects. 

2. Validate the rotor response calculation with measured data obtained from a model rotor tested 
in a wind tunnel. 

3. Formulate the higher harmonic control analysis for open-loop and closed-loop controllers. De- 
termine the power requirement for hydraulic servo-actuators. 

4. Perform systematic parametric studies to investigate the effects of blade design parameters on 
the vibratory hub loads, HHC inputs, and power requirements. Results from these studies might be 
used to assess the role of blade dynamics on the efficiency of HHC systems and explain the differences 
in HHC amplitude requirement between different rotors. The design parameters include blade stiffness 
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distributions in the flap, lead-lag, and torsion directions, chordwise offset of the blade center of mass 
from the elastic axis, chordwise offset of the elastic axis from the blade quarter chord, and blade linear 
pretwist. 

5. Examine application of HHC for articulated and hingeless rotors, in particular the HHC ampli- 
tude and actuator power requirements. 

6. Investigate effects of HHC on blade bending stresses, control loads, and rotor performances 
under different flight conditions. 

7. Investigate the application of HHC at severe flight conditions including high speeds and high 
thrust flight conditions, in particular the influence of higher harmonic blade pitches on blade stall. 

The analytical study is conducted using a comprehensive coupled aeroelastic analysis developed 
at the Rotorcraft Center, University of Maryland (refs. 77-80). The formulation is based on the finite 
element method in both space and time. In this formulation, the rotor blade is modeled as an elastic 
beam undergoing flap bending, lead-lag bending, torsion, and axial deflections. The spatial finite element 
method based on the Hamilton principle is used to discretize the blade into a number of beam elements, 
each having 15 degrees of freedom. 

Since accurate prediction of the blade airloads is crucial in the determination of the aeroelastic 
responses and the higher harmonic responses, a time domain nonlinear unsteady aerodynamic model 
developed by Leishman and Beddoes (refs. 81 and 82) has been used. This model consists of an 
attached compressible flow formulation (linear model) along with a representation of the nonlinear 
effects due to flow separation and dynamic stall. The impulsive — or noncirculatory — loads and reverse 
flow corrections are also included. Calculation of the rotor nonuniform induced inflow is performed 
using a free wake analysis adapted from a Comprehensive Analytical Model of Rotor Aerodynamics 
and Dynamics (CAMRAD) (ref. 83). The blade nonlinear periodic responses are solved from the 
normal mode equations using a temporal finite element technique. The vehicle trim controls and the 
blade responses are calculated iteratively as one coupled solution using a modified Newton method. 
Calculation of the HHC responses forms an integral part of the steady response solution in which the 
basic trim controls are superimposed with the higher harmonic components. 

The HHC algorithms are formulated in the frequency domain. A linear quasi-static transfer matrix 
is adopted to relate the harmonics of HHC inputs to the harmonics of vibratory hub loads, and both 
local and global models are considered. The optimal HHCs are calculated based on the minimization 
of the vibratory hub loads expressed in term of a quadratic performance index. The transfer matrix 
is computed off line, and, for the local model, it is updated based on the secant method. Using the 
multiblade coordinate transformation, the HHC inputs in the rotating frame are transformed to the fixed 
frame, in terms of the swashplate higher harmonic motions. The higher harmonics swashplate motions 
are generated using a number of servo-actuators, and the hydraulic power required by these actuators 
is then calculated for a specified HHC input. 
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2 FORMULATION 


This section is divided into two parts: the first part describes the baseline analysis with no HHC 
effects, and the second part describes the formulation of the HHC analysis. 

In the first part, the coordinate systems and the Hamiltonian formulation of the blade governing 
equations are introduced. Then the blade aerodynamic and the rotor induced inflow models used in 
the analysis are described. The next section shows the calculations of the blade and rotor hub loads, 
and this is followed by the description of the coupled trim analysis. The subsequent section presents 
the spatial finite element method used in the formulation of the blade governing equations. Then, the 
solution technique for the blade responses as coupled with the vehicle trim analysis is shown. 

The second part, starting with section 2.9, describes the higher harmonic control analysis. First, 
the HHC model and the control algorithms used to suppress the rotor vibratory hub loads are presented. 
These are followed by the formulation of the actuator power requirement for the HHC system. Then, 
the implementation of HHC into the overall analysis is discussed. 


2.1 Coordinate Systems 


The helicopter model consists of a rigid fuselage and a main rotor having N b blades. The Cartesian 
coordinate system (A, Y, Z) is attached to the rotor hub which remains fixed to the fuselage, as shown 
in figure 2.1. This coordinate system, with the origin at the rotor hub center, represents the fixed or 
nonrotating frame. The Z axis is along the shaft axis, and the X and Y axes are in the hub plane, 
perpendicular to the shaft axis. 

The rotor blade is modeled as an elastic beam rotating at a constant angular velocity Cl: Figure 2.2 
shows the blade coordinate system. A Cartesian coordinate system (x, y, z) is attached to the undeformed 
rotating blade, which is at a precone angle /3 p with respect to the hub plane. This coordinate system, 
with the origin coinciding with that for the fixed frame, represents the rotating frame. The x axis 
coincides with the blade elastic axis and points outboard; the y axis is in the plane of rotation and points 
toward the blade leading edge; the 2 axis is tilted at an angle (3 p relative to the shaft axis. 

For elastic deformation, a point P on the blade elastic axis undergoes displacements u, v, w in 
the x, y, z directions, respectively, and occupies a position P f on the deformed elastic axis. As shown 
in figure 2.2, the cross section containing P' undergoes a rotation 9\ about the deformed elastic axis, 
where 6\ can be expressed as (ref. 84) 

01 = 0 + 0 ( 2 . 1 ) 

and the twist 0 is given by 



where 0 is the elastic twist of the blade cross section about the deformed elastic axis, 0 is the geometric 
twist, 6 is the control pitch angle, v' = dv/dx is the lag slope, and w" = d 2 w/dx 2 is the flap curvature. 
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An orthogonal coordinate system (f , 77, Q is introduced at the deformed position P' such that the 
£ axis is tangential to the deflected elastic axis. The 77 and ( axes are the principal axes of the deformed 
cross section and remain normal to the deformed elastic axis. This deformed blade coordinate system 
is also shown in figure 2.2. The assumption, which is based on the fact that the blade cross section 
remains normal to the elastic axis after deformation, is in accordance with the Bernoulli-Euler beam 
theory. In general, the plane containing ( 77,0 and that containing (y,z) are not parallel. 

2.2 Hamiltonian Formulation 


The variational method of the Hamilton principle provides a way for finding the equations of 
motion of a dynamical system, and such approach has been adopted in the analysis. For a conservative 
system between prescribed initial and final conditions at time t\ and £ 2 * the Hamilton principle states 
that the actual motion of a dynamical system is that particular motion for which the time integral of the 
difference between the potential and kinetic energies is minimum (ref. 85). Since aeroelastic systems 
are nonconservative, the Hamilton principle is generalized to include forces which are not derivable 
from a potential function. The resulting generalized Hamilton principle, applicable to nonconservative 
systems, can be expressed in variational form as 



-ST- 6W) dt = 0 


(2.3) 


where ST and SU are, respectively, the variations of kinetic and strain energies, and SW is the virtual 
work done by external forces. 


From reference 84, the variation in strain energy for a Bernoulli-Euler beam undergoing moderate 
deflections in the axial, lead-lag, flap, and torsion directions is 

SU = (F(Su 4- v'Sv' + w'Sw') 

+ (GJ<f> + EAk\(e + <t>)'(u' + \v n + i^' 2 ) 

Z Z 

+ EBxO?(f)' - EB 2 9[(v" cos 6 X + w" sin 9^)8 ft 
-I- (EC\(j)" -f EC 2 (w" cos 9\ - v" sin 0\)8(j> f 

+ ((M z - EAe a {u + \-v’ 2 + i«/ 2 ) - EB 2 0'(f) f ) cos6\ 

+ (My - EC\(j) n ) sin 0i){8v" + w"6lj>) 

+ {{M z - EAe a (u! + ^v' 2 -f ^ /2 ) - EB 2 Q'$) sin#! 

— (My — EC\4> n ) cos 9i)(Sw ,f — v"S4>)) dx (2.4) 

where the axial force F is defined by 

F — EA(u ! + -v’ 2 + -w' 2 + — e a (v" cos 9\ -f w" sin #i)) (2.5) 

Z Z 
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and 


( 2 . 6 ) 

(2.7) 


My = EI y (v" sin 0\ — w" cos 0\) 

M z — El z (v" cos 6\ +tL> // sin^x) 

and the sectional constants are 

a= JJ a dvdc 

Ae a = fj r)dr)d( 

7 »= JJa (2dvdc 

Iz= JL n2drid( 

a *&= JJ A (v 2 +e)d V d< 

J = JJa + + 

B l = f J A iv 2 + ( 2 ) 2 dr} d( 

B 2 — JJa + C 2 )^^C 

°i = JJa x2d7]d ( 

°2 = fj A (><dr)d( 

where A is the blade cross-sectional area effective in carrying tension; e a is the offset of the tensile 
axis from the elastic axis, positive forward; Iy and I z are the flapwise and chordwise area moments of 
inertia of the cross section; Ak ^ is the cross-sectional polar moment of inertia; J is the torsional constant 
including cross-sectional warping; A is the warping function which is assumed to be antisymmetric about 
the chord line; C x is the warping rigidity; and B x , B 2 , and C 2 are other sectional constants. Furthermore, 
it is assumed in the analysis that the structural properties of the cross section are symmetric about the 
chord line. 
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The variation in kinetic energy for a rotating elastic beam is 

fR 0 

ST — J m((frz + 2Qf> - u)6u 

+ (ft 2 (v -f- cos^i) + 2 Q((3pw — u) 

+ 2e 5 £l(r/cos#i + w'sin Q\) —v + e^sin 
(/?p(£2 2 x -J- 2Qu) + ti) + e g (j> cos )6iy 

“ ^ (^m 2 — ) cos ^isin 6 \ 

+ e fl fi 2 x(u/cos0i - t/sin 0\) + egtfvsinOi 
+ e g Q 2 f3 p x cos 0\ — e g (vsinOi — wcos0i))6$ 

- e g (Vt 2 x cos $i + 2 ttv cos S\ Sv' 

— e^(f2 2 j;sin 6\ + 2fTri;sin 6\ )6w') dx (2.8) 

where the sectional integrals are defined by 

m = JJ A Psdr]d( 
me g = JJ A PsP dr/ d( 
m ^i = JJ A Ps ( 2 dr] d( 

mk m 2 = JJ A PsP 2 dr] d( 

1.2 _ r 2 i t.2 

— ^777,2 

In the above expressions, m is the blade mass per unit length; p s is the blade structural density; eg is 
the center-of-mass offset from the elastic axis, positive forward; mk 2 ll and mk are the flapwise and 
chordwise mass moments of inertia per unit length; and k m is the cross-section radius of gyration per 
unit length. 


An ordering scheme has been used in the derivation of the variations in potential and kinetic 
energies. Basically, the ordering scheme provides a set of guidelines to neglect the higher order terms 
which appear in the nonlinear equations. These guidelines are employed to simplify the equations of 
motion. By introducing e, with e <C 1, a parameter of order v/R or w/R, the assumed orders of 
magnitude for the nondimensional quantities are 


EA/m 0 ft 2 R 2 0(e~ 2 ) 

x, m/m 0 , 0 . 75 , 6\ c , 6\ s 0 ( 1 ) 

EIy/m 0 ft 2 R 4 , EI z /m 0 Q 2 R 4 , GJ/m 0 ^ 2 R 4 0(1) 

v/R, w / R , (f>, (5p, kj±/ R , k m\ / R, k m 2 /R, EB2/rn 0 td 2 R 5 (9(e) 

e g /R,e a /R,e d /R 0(e 3 / 2 ) 

u/ R , EB \ / m 0 Q 2 R Q 0(e 2 ) 
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In general, terms with the highest order for a particular equation are considered first order terms. 
First order terms are retained in the u equation. For instance, if the highest order terms in the u 
equation have an order of then all terms with order are retained, while higher order terms such 
as £ 3 ,£ 4 ,... are neglected. Similarly, second order terms are retained in the v and the w equations, 
while some third order terms are retained in the $ equation. Furthermore, these guidelines should 
be applied within the total energy expression, and the terms retained should constitute the self-adjoint 
characteristic of the structural matrices. That is, the resultant blade structural mass and stiffness matrices 
should be symmetric, while the antisymmetric structures of the damping matrix should reflect gyroscopic 
couplings. 

The expression for 6W is 

- J (L^6u + LySv + L^bw + M^<50) dx (2.9) 

where and are the distributed airloads in the x,y,z directions, respectively, and is 

the sectional aerodynamic pitching moment about the deformed elastic axis. Based on the Hamilton 
formulation, a finite element analysis is formulated in section 2.7. 

23 Aerodynamic Loads 


This section describes the calculation of the blade section airloads. The external loads L £ 

and in equation 2.9 are the distributed airloads acting on a rotor blade. Calculation of the blade 
airloads requires information about the relative flow velocity and the Mach number as seen near a blade 
section. Therefore, the derivation of the relative flow velocity is presented first, and this is followed by 
the presentation of the aerodynamic model. 

2.3.1 Derivation of Blade Velocity 


In particular, the resultant velocity seen at a blade section consists of the incoming flow velocity, 
the instantaneous blade motions, and the induced inflow due to the trailed and shed vorticities. The 
general expression for the resultant blade velocity at a radial station x is given by 

V = -V w + V 6 (2.10) 

where \ w is the wind velocity with contributions from the vehicle forward speed and the rotor inflow, 
and V 6 is the blade velocity relative to the hub fixed frame resulting from blade rotation and blade 
motions. The expression for the wind velocity is 

V™ = (fiUR)I - ( XUR)K (2.11) 

where fi = V cos a s /QR is the rotor advance ratio; V is the vehicle forward speed; a s is the rotor 
shaft tilt, positive forward; A is the rotor nondimensional inflow; and UR is the rotor tip speed. The 
unit vectors (/, J,K) are the basis vectors for the nonrotating coordinate system. The rotor inflow A 
consists of two components and is expressed as 

A = fi tan a s + A; (2. 1 2) 


23 



where /i tan a s (fig. 2.3) is the component of the forward flight velocity perpendicular to the hub plane, 
and Xi is the nondimensional rotor induced inflow associated with the lift on the rotor. Calculation of the 
induced inflow A 7 involves extensive formulation, and the details are deferred to the next section. Since 
the velocity components in equation 2.11 are expressed in the nonrotating system, two transformations 
are required to change these to the rotating undeformed frame. The first transformation, which converts 
the coordinates in the hub fixed frame to the rotating frame with no precone, is given by 

cos’J' sinT' 0 " 

T 1 = -sin\£ f cos^ 0 (2.13) 

0 0 1 J 

The second transformation relates the rotating coordinate system with no precone to that at a precone 
angle j3 p and is expressed as 

cos Pp 0 sin (3 p 

t2 = o 10 (2.14) 

— sin (3 p 0 cos (3 p 

Therefore, subsequent substitution of and in equation 2.11 results in the wind velocity expression 


Vw — V X l 4- Vyj + V z k 

(2.15) 

where V x ,V y , and V z are the components of the wind velocity, and 

V x ~ fiQR cos # cos (3 p - Af }R sin (3 p 

(2.16) 

V y = — /iO R sin 

(2.17) 

Vz = — fiQR cos — Af 2 i?cos (3 p 

(2.18) 

The general expression for the blade velocity is 

Vpr + flxr 

(2.19) 


where r — [x\, y\, 2 q), r — (x,y,z), and Q — (ST sin 0, Qcosftp). In particular, the expressions for 
r and r consist of the blade motion terms, and, being calculated at an arbitrary position ( 77 , £) of the 
blade cross section, \ b takes the form 


Y b = VbJ + VbyJ + V hz k ( 2 . 20 ) 

where V bx ,V by , and V bz are the components of the blade velocity, and 

•j 

V bx = u — \<f) — (7/ + w' 6\)(r) cos $i — £sin< 9 i) 

- (w' - v'Gi) (//sin 61 + £ cos 0 \ ) — Q(v + 77 cos 0\ ~Csin<9i)cos (3 p (2.21) 

V by = v - $i (77 sin 61 + ( cos 9 \ ) 

+ Q((x + u — \<j)' — 7/(77 cos #1 ~ 

- 77/(77 sin 9 \ + £ cos 9 \ )) cos /3 p — (w 4 - ysin 9 \ + 77 cos 9 \ ) sin ( 3 p ) (2.22) 
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(2.23) 


V bz — w + # 1(77 cos 0i - Csin 0 i) + Q(v + r?cos0i - Csin0i) sin fi p 

The resultant velocity at the blade cross section, expressed in the rotating undeformed coordinate 
system, is 

Y = U x i + U y j + U z k 

= (V bx - V x )i Hh (V^ - V y )j + (V bz - V z )k (2.24) 

Since a strip analysis is used, airloads (the aerodynamic loads) are computed based on the angle of 
attack at the three-quarter chord location. With respect to the blade coordinate system, this requires the 
velocity components to be calculated at (0, r} r ). As such, the resultant velocity expressions simplify to 

U x = u — (v' -f w'0i)r] r cos 0 i — ( w ’ — i/ 0 *i)? 7 rsin 0 i 

— 0 (u + ? j r cos 0i ) — /j,£IR cos 4/ + Af }RP p (2.25) 

U y — v — 0i ? 7 r sin 0i + [iQR sin 

+ Q(x + u — vr} r cos 0i - w'rjrsin 0i — (w + r ) r sin 0 i )/3 p ) (2.26) 

U z — w + 0*i 77r cos 0i + tt(3 p (v + 7] r cos 0i ) + fiQR/3 p cos *1/ + A FIR (2.27) 


In deriving the above expressions, small angle assumptions were made on the precone angle 
that is, sinfip w f3 p , and cos (3p « 1. 

Finally, the blade velocity components are transformed into the blade deformed frame. Such 
information is required to calculate for the blade section airloads in the deformed frame. In particular, 
the transformation relating the coordinate system from the undeformed to the deformed frame is given 
by (ref. 84) 

v 1 w' 

(1 — tj-) cos 0i — v'w'sin 6\ sin0i(l — 

/2 /2 
-(1 — -^-)sin0i — r/u/cos0i cos0i(l — ^ 2 “) 

Thus, the resultant blade velocity expressed in the rotating deformed frame is 

V = Upe £ + Upe p + Upe £ (2.28) 


rji3 _ 


/2 

V 

~T 


w f2 

~Y~ 


— ( v f cos 0i -f u/sin 0i ) 
u'sin 0i — w' cos 0i 
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where 


— = u — v + v'(x + /x sin '!>) - /xcos^(l - (3 p w f ) + \(/3 p + w’) 

— r) r cos 9(1 + v f — 9w') + r] r sin 6 ((f) — w' — 9v') 

4- v'v -f w'w + -ficos fy(v' 2 + w' 2 ) 

& 


Ur 

SIR 


Up 

SIR 


cos 9((v + w) — i<;/3p + 0(A + w) 4- r/(u — w) 
v' 2 

+ (x + /xsin\I>)(l - — ) -f/xcos^v' + 0(/3 p + w f ))) 

+ sin <9((w + A) + v(/3 p + w) - <f>(v + u) - iiw 

- (x 4- (is in ^)(v'w' + + M cos ^(w' + fi p — 4>v')) 

cos#((?i; + A) -f [3 p v -f (v — u)w f + fxcos^(p p + w/ — </>?/) 

- (x + /x sin + </>)) + sin <9(-(v + w) + (tx - u)i/ + iu/? p 

/ - t / 2 

- (f)(w + A) - fi cos $f(v + (f)((3 p + w)) - (x + fi sin #)(1 )) 

2 

4- r} r (9 -F 4> + w' + (3 P + /xcos #0*^) 


(2.29) 


(2.30) 


(2.31) 


In the above expressions, only the second order terms are retained, while the higher order terms 
are neglected. This assumption is consistent with the fact that the blade deflections and velocities are 
small compared to other terms. Furthermore, small angle assumptions were made on <£ such that 

sin(0 + (j>) « sin6> + 4>cos9 (2.32) 

cos(0 + (j>) « cos 9 — sin 9 (2.33) 


2.3.2 Aerodynamic Modeling 


Calculation of the blade section airloads is performed in the time domain using a nonlinear unsteady 
aerodynamic model developed by Leishman and Beddoes (refs. 81 and 82). Elements of the model and 
its application to the aeroelastic analysis are described in this section. Basically, this model consists of 
an attached flow formulation along with a representation of nonlinear effects due to flow separation and 
dynamic stall. 

The basic approach adopted in formulating the unsteady aerodynamic model is to identify the 
physical phenomena associated with each process and to model the whole process by using a combination 
of linear dynamical systems. Basic elements of the model, representing each process, are described in 
terms of the indicial response functions. For a linear system, the indicial response is the system 
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response to a unit step input, with zero initial conditions. If the indicial response is known, then the 
system steady response can be computed using the Duhamel integral. In particular, if <j)(s) is an indicial 
response function characterized in terms of the time variable s, then the system response y(s) to an 
arbitrary forcing function u(s ) is 

y(s) = <t>(s)u{ 0) +/</>(«- a )~^- da (2.34) 

VO acr 

where a is the dummy time variable of integration. Since the rotor aeroelastic analysis is primarily 
concerned with the periodic force response, the transient response associated with the initial condition 
is neglected. 


Attached flow formulation- The attached flow formulation is applicable to flow regime where the 
airfoil lift and pitching moment vary linearly with the angle of attack. Even though nonlinear effects 
associated with the stall behavior have not yet been included, this formulation is a prerequisite to the 
complete model including flow separation and dynamic stall. The attached flow formulation has been 
generalized for the compressible flow, and, furthermore, it has been validated for airfoils that undergo 
arbitrary motions up to the transonic flow regime (ref. 86). 


For the range of angle of attack typically encountered in rotor analysis, the following assumptions 
are made (fig. 2.4) 

L w « L (2.35) 

with 

C L = C 0 + Ci a (2.36) 

where C\ is the compressible lift curve slope, with value taken from the static airfoil data at M = 0.3, 
and scaled by the Prandtl-Glauert correction factor giving 

c, = ^ ] f d™. (2.37) 

where /? = Vl - M 2 . Furthermore, since the angle of attack a is small 


“ = tan_1 {-w 

-Up 

~ -jf- (2.38) 

where Up is the downwash at the three-quarter chord location, and Up is the incoming velocity com- 
ponent, parallel to the chord line. 


With these assumptions, L w can be expressed as 

Zw(s) = pcUfeo + L c w (s) + Li(s) (2.39) 

where L%j is the circulatory lift, and L l w is the noncirculatory lift. In particular, T w is the component 
of the airloads acting perpendicular to the chord line and is usually referred to as the normal force, 


27 



also denoted by N. Using the indicia! function formulation, the component of circulatory lift can be 
expressed as 

L w ( s ) = ~ 2 P cC iUt <t>c(s)U p(ti) + (f) c (s - do- . (2.40) 

The circulatory indicial response (j> c (s) is given by (ref. 81) 

Ms) = 1 - [Aiexp(-6 1/ S 2 s) + J 4 2 exp(-6 2 /3 2 s)] (2.41) 

where the terms in bracket represent the lift deficiency due to vortices shed downstream into the 
airfoil wake. Without the Mach scaling parameter (3 2 , this indicial function is similar to the classical 
Kussner function and differs only through the values of the coefficient. For this indicial function, 
A\ = 0.3; A 2 = 0.7; b\ = 0.14; 62 = 0.53, while the coefficients for the Kussner function are (refs. 60 
and 61): A\— A 2 = 0.5; b\ = 0.13; 62 = 1.0. The variable s represents the distance an airfoil travels 
in terms of semi-chord, and it can be expressed in terms of rotor coordinate as 

• s = £ r utw < 2 - 42 > 


Similar to the circulatory lift, the impulsive lift component is given by 


'( 5 ) = ~2 PCU T ^(^^(0 ) + f Q & i( s ~ a ) 


dU P (o) 


(2-43) 


where <fo(s) is the impulsive indicial response function. This function is derived based on the dynamics 
of the noncirculatory loading due to pressure wave propagation, and the initial value is computed using 
piston theory (ref. 61); that is 

Ms) = T7 exp (2.44) 


M s ) — exp \2MK^J 

where K a is a time constant dependent on the Mach number (ref. 82), and 




2(1 - M) + 2irpM 2 {A l bi + Mb 2 ) 


(2.45) 


This impulsive loading can be considered as the compressible analog of the apparent mass terms used 
in many incompressible analyses. 

Following the formulation of reference (ref. 87), the chord force L v , also denoted by C, can be 
expressed as 


— 1 2 

L v — rj a L sin a - -pUpcdo cos a 


(2.46) 


where % is the leading suction recovery factor which accounts for the inability of the airfoil to realize 
100 percent suction in a real flow, and do the viscous drag or skin friction drag component. 

With regard to the more familiar drag coefficient expression, which is often presented in terms of 
a quadratic polynomial in a as 

Cd = do + d\a. + d2oP‘ (2.47) 
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it is worth noting that the drag coefficient can also be expressed in terms of the lift and chord forces; 
that is 

Cd = do + C L sin a ~C c cos a (2.48) 

with C c being the chord force coefficient. It can be shown that these two descriptions are equivalent, 
with the latter being more suitable for unsteady calculations. Assuming that the a is small, then 
equation 2.48 can be approximated as 


Cd ~ do + (1 - r] a )C L a (2.49) 

where the second term in the above expression is actually the pressure drag component. Together with 
equations 2.36 and 2.47, it can be shown that 

d\ = (1 - v)Co 

d>2 = (1 - v)Ci 

These results imply that the pressure drag can be expressed solely in terms of the normal force. As 
such, the unsteady pressure drag can be computed from the time-dependent lift and chord forces. 

The radial drag L u is computed based on the fact that only the viscous drag is affected by the 
velocity component in the radial direction (ref. 88). Hence, the pressure drag does not contribute to the 
radial drag component. In the radial direction of the rotor blade 

Lu — Dr — 2 p{Ur + Uj) c do sin A 
where A is the blade yaw angle, and 

• a - u r 
sin A = . = 

M 

With the above expression, L u simplifies to 

L u = ^pcd 0 Jul + U$,U R (2.50) 


As in the case of lift, the unsteady aerodynamic pitching moment also consists of circulatory and 
impulsive components. However, the pitching moment is dependent on the modes of forcing, or the 
types of airfoil motion. These modes of forcing basically govern the types of pressure distribution on 
the airfoil, which, in turn, have a direct effect on the pitching moment. In the analysis, the airfoil angle 
of incidence is assumed to consist of two types of airfoil motions: a pure heaving or vertical motion 
along the £ axis, denoted by h, and a pure pitch rate or the airfoil pitching motion about the elastic 
axis, denoted by q . Together these modes of airfoil motion contribute to a as 


where 


a 


j}_ q 

U T 2 


—2rj r a 

~U T 
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and 7 j r is the offset of the airfoil pitch axis from the three-quarter chord location. 

With these modes of forcing, the unsteady pitching moment about the airfoil elastic axis is given 

by 

M <t> = \puh 2 fa ~ e d cL c w + M' h + Mtj + M' q (2.5 1 ) 

where f 0 is the pitching moment coefficient at zero lift about the aerodynamic center, and e d is the 
chordwise offset of the aerodynamic center behind the elastic axis, expressed in terms of chord. The sec- 
ond term in the above equation represents the component of the circulatory lift acting at the aerodynamic 
center. The third term represents the contribution of the airfoil heaving motion to the noncirculatory 
pitching moment and is expressed as 


M h ~ ~2 pc 


4>M h (s)h(0) + J o 4> l Mh ( s -ex) ^ 

with the corresponding indicial response function given by 


dh(cr) 


do\ 


(2.52) 


^M h ( 5 ) - ^3 exp ( — - , S — ) + A 4 exp ( o S — 
V 2M&3 Km h ) \ 2M64 Km h 

where A 3 = 1.5, A 4 = -0.5 ,b x = 0.25 ,b 2 = 0.1, and 

A364 + ^463 


K Mh = 0.8 


6364(1 - M) 


(2.53) 


Physically, this indicial function represents the decay of the impulsive loading along with the forward 
movement of the center of pressure. Due to a unit step input in A, the pressure loading is initially 
uniform across the chord length, and the center of pressure is at midchord. As the loading decays 
exponentially, the pressure distribution becomes asymmetric, and the center of pressure tends toward 
the airfoil quarter chord. 


The fourth term in equation 2.5 1 is the circulatory pitching moment due to the airfoil pitch rate q 
about the elastic axis and is given as 


A/£(s) = 


7T 


1 6/3 


pc 2 U T 


0M,W<Z(O) + f 4> c Mq {s — a )~r~ d/r 


(2.54) 


This is usually referred to as the pitch rate induced camber term, and its quasi-steady value is taken 
from the classical thin airfoil results. The associated indicial response is 


= 1 - exp (— &5/3 2 s) 

where 65 = 0.5. 


(2.55) 


The noncirculatory pitching moment, the last term in equation 2.51, has contribution from the 
airfoil pitch rate q and is expressed as 


Ml(*) = -~P<?U T 


<f>M q (s)q(0) + J <l? l Mq (s - da 


(2.56) 
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The corresponding indicial response given by 


%,( s ) 12M exp \ 

where 

K 5.6 

M<? 15(1 - M) + 3tt/?M 2 & 5 


(2.57) 


Separated flow formulation- Extension of the aerodynamic formulation to the nonlinear regime 
requires the modeling of two distinct phenomena associated with the airfoil stall behavior, namely 
trailing edge and leading edge separations. The phenomenon of progressive trailing edge flow separation 
is associated with most types of airfoil stall under both static and dynamic conditions. Trailing edge 
separation starts from the airfoil trailing edge, and with increasing angle of attack, the separation point 
progresses toward the leading edge region. Such effects contribute to the nonlinear behavior in lift, 
drag, and pitching moment due to the loss in circulation. If supercritical flow develops under dynamic 
conditions, then the stall process is initiated by the leading edge or shock-induced separation. This 
phenomenon is characterized by the shedding of a strong vortex from the leading edge region. The 
leading edge vortex induces a large pressure disturbance wave, moving aft on the airfoil upper surface; 
produces abrupt changes in the flow field; and contributes to the large transient lift, drag, and pitching 
moment. This is a catastrophic process known as dynamic stall. 

To simulate the nonlinear effects associated with trailing edge separation, Kirchhoff’s flow model 
is used (ref. 89). This model relates the airfoil static normal force (or lift) and chord force behavior to 
the trailing edge separation location, a location at which the flow starts to detach from the airfoil upper 
surface. In this model, the expressions for the normal and chord forces are 

Cn = Cl a (2.58) 

C c = riaCi^a 2 (2.59) 

where / is defined as the chordwise distance between the flow separation point and the airfoil leading 
edge, expressed in terms of chord and shown in figure 2.4. In particular, the value of / varies between 
zero and one and dictates whether the flow is fully attached (/ = 1) or completely separated (/ = 0). 

A general expression for the pitching moment was not provided by the Kirchhoff model, and, 
therefore, the following empirical relation is used (ref. 81) 


Cm = /o + [&o + &i (1 - /) + &2 sin Cn 


(2.60) 


where ko — \ — x ac is the aerodynamic center offset from the airfoil quarter chord, and k\ and &2 are 
empirical constants obtained from a curve fit to the static pitching moment coefficient for a particular 
airfoil. 
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Implementation of the above relations (eqs. 2.58 and 2.59) requires the explicit representation of 
the airfoil separation point / in terms of the angle of attack. In particular, using equation 2.58, an 
effective separation point can be inferred from the static lift coefficient for a given airfoil at a given 
Mach number. Specifically, if the static lift coefficient is known, then 


/ = 



(2.61) 


where C\ is the lift curve slope at zero lift value. Curve-fitting relations are then used to describe the 
variation of / with a, and as such 


fl 0.3 exp ^ for a < a\ 

| 0.04 -f 0.66 exp ^ for a > a\ 

where a\ is the angle of attack at which / = 0.7, a breakpoint dividing flow behavior into partially 
attached (f > 0.7) and fully separated (f < 0.7) conditions. For a pair (a,/*) and oq, the empirical 
constants Si and S 2 can be calculated by inverting equation 2.62 as 


1 

si 


1 

S~2 


1 

a — aq 
1 


In 


aq — a 


In 




l 0.3 j 

(2.63) 

/ - 0.04 \ 


” 0.66 ) 

(2-64) 


For a given Mach number, the variation of / with a is characterized in terms of only three empirical 
parameters: aq, Si, and S' 2 . Thus, this model provides a compact structure for the reconstruction of 
the airfoil static nonlinear behavior. If the effects of the airfoil unsteady pressure distribution and the 
onset of separation can be evaluated, then this model can also be extended to the dynamic regime. 

To model the onset of separation, it has been established by Beddoes (ref. 81) that only the peak 
pressure at the leading edge or ahead of the shock wave is important. Therefore, generalization of 
the dynamic effects based on the critical aspect of pressure distribution provides a useful means of 
implementing the criteria for leading edge or shock-induced separation. This generalized criterion 
also provides a way for finding the unsteady trailing edge separation based on the assumption that the 
boundary layer and separation development are also dictated by the behavior of the leading edge pressure. 
Hence, implementation of the criterion would modify the trailing edge calculation and simultaneously 
initiate the leading edge vortex shedding. 

Under unsteady conditions, position of the separation point is governed by the dynamics of both 
the airfoil pressure distribution and the boundary layer responses. First, it is required to relate the airfoil 
unsteady leading edge pressure response to the forcing. Based on unsteady airfoil data, the dynamics 
of the airfoil pressure distribution with respect to the lift forcing can be represented by a first order 
linear dynamical system with a time constant T p . In particular, the system input is C^ L , the airfoil 
attached flow normal force coefficient, and the system output is C' N , the substitute normal force value. 
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The indicial response of C' N is given by 


<Pp{s) = 1 - exp 



(2.65) 


where 

^ *<*> 

and L w is given by equation 2.39. The time constant T p depends on the Mach number and can be 
determined from unsteady airfoil data. 


For the unsteady case, the response of the boundary layer is also subjected to dynamic behavior. 
Again, modeling of this process is accomplished by representing the boundary layer response by an- 
other first order linear dynamical system with a Mach dependent time constant 7y. In the model, the 
quasi-steady separation point, denoted by /', is computed using the substitute normal value. This is 
accomplished through an effective angle of incidence ay given by 


The resulting value of /' is evaluated by substituting the value of ay into equation 2.62. Then the 

indicial response of / , the effective trailing edge separation point, to changes in the airfoil pressure 
distribution takes the form 


<pf{s) = 1 — exp 



(2.67) 


In particular, / represents the pressure distribution response, while f ,f includes additional effects of the 
boundary layer response. Both values of Cjy and f" can be computed using the superposition integral 
in the form of the Duhamel integral, and the effective separation point is then used to calculate the 
airfoil nonlinear lift, drag, and pitching moment. The value of the time constant Ty can be obtained 
from unsteady boundary layer calculations. 


Physically, dynamic stall is characterized by the leading edge separation and the shedding of 
concentrated vorticity from the leading edge region. This vorticity is swept downstream over the chord 
and produces a significant increment in the airfoil lift, drag, and pitching moment. Since dynamic stall 
is a catastrophic process, the stall onset is modeled using a criterion such that the separation is initiated 
when the substitute normal force C' N attains a critical value, C Nv Based on the calculation of C' N 
proposed earlier, this criterion implies that, under dynamic condition, the critical lift condition may be 
reached at a higher angle of attack than the static condition. Typically known as dynamic stall delay, 
this effect contributes to the delay in the occurrence of stall under dynamic conditions. With regard 
to the angle of attack forcing, the dynamic stall delay is represented by the unsteady effects associated 
with the attached flow normal force and the time constant T p . Also, it is worth pointing out that the 
value of T p depends on the Mach number. The value of C corresponds to a critical leading edge 
pressure and is closely related to the airfoil maximum static lift coefficient. Furthermore, the value of 
C Nl depends on both the airfoil geometry and Mach number. 
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Once leading edge separation has been initiated, the vortex lift C NV , being the excess lift induced 
by the dynamic stall process, is governed by the incremental lift quantity Cy. The value of C v is 
related to the difference in lift between the attached and the separated flow regimes. In particular, this 
process is governed by the rate of change of Cy. In the limit as the rate of change tends to zero, the 
value of Cjyy decays exponentially to zero. Modeling of this process is realized by yet another linear 
dynamical system of the first order with a time constant Ty. Employing this model, the response of 
Cny to the time rate of change of Cy can be expressed using the convolution integral as 

C N v(s) - <t>v{s)Cy( 0) + [ (j> v (s - cr) ~¥~^ da (2.68) 

./0 da 

where the indicial response function <f> v (s) has the form 


and 


(fiv(s) = exp 



Cy = C Nl 1 



(2.69) 


(2.70) 


As the leading edge vortex is swept downstream, the airfoil center of pressure varies significantly 
with the chordwise location of the vortex. The vortex convection rate is monitored using a nondimen- 
sional time based on the distance traveled by the airfoil in semi-chords. This rate has been determined 
from experimental data to be roughly less than half the free-stream velocity, while the actual rate is 
dependent on the Mach number. Based on physical considerations, a suitable representation of the airfoil 
center of pressure behavior takes the form 

c/v = i(i- cos (iz)) 

where r v is the vortex time, and Tyi represents the vortex convection rate. At the onset of separation, 
the vortex is at the leading edge, and r v is equal to zero. When the vortex reaches the airfoil trailing 
edge, r v is equal to Tyi. Therefore, by monitoring the excess lift and the center of pressure location, 
the vortex induced pitching moment can be calculated by 


Cmv — ~CPy C^y (2.72) 

The values of Ty and Tyi are determined from the airfoil unsteady data. After the vortex passes 
over the trailing edge, the effects of the vortex lift on the airfoil decay rapidly, and this phenomenon is 
accomplished by halving the time constant Ty for the period Tyi < r v < 2 Tyi. 

Further modifications to the basic separated flow formulation have been made to better capture the 
interactions between the dynamic stall process, the trailing edge separation behaviors, and the airfoil 
kinematics. In most cases, these interactional effects are represented by modifying the appropriate time 
constants associated with each process. Specifically, increasing a time constant has the effect of slowing 
down the process, whereas reducing the time constant produces the opposite effect. Further details on 
these modifications can be found in reference 81. 
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A final modification made to the basic model concerns the remodeling of the Kirchhoff chord force 
expression (i.e., eq. 2.59). Such modification is necessary to model the drag divergence induced by the 
flow behavior after the onset of gross separation. For this, an additional parameter has been introduced 
to the chord force expression, which becomes 

Cc = VaCiyjf'a 2 ® 

where 

$ = ~ C Ni) 

for C' N > Cjsf l . The actual value of Kj can be determined from the airfoil static drag data. 

In forward flight, the rotor is also subjected to a reverse flow region that increases with increasing 
forward speed. The reverse flow region is characterized by a complete reversal of the inplane velocity 
vector as shown in figure 2.5. This effect contributes to the reversal of the chord force and the pitching 
moment, and causes the shifting of the aerodynamic center from the region near quarter chord to 
approximately the three-quarter chord location. The lift and radial drag components are independent of 
the direction of the inplane velocity vector. 

In the analysis, the computation of the blade airloads is performed using a finite-difference ap- 
proximation to the Duhamel integrals. And with the inclusion of reverse flow correction, the unsteady 
aerodynamic loads, at the rotor time sample ^ n , are expressed in a discrete-time form as 

Lu('S’n) = ~y [dosJU 2 R + U%.U R 

Lv(*n) = -ysign (U T ) d 0 U f. - 'naC 0 \ff"$\U T \Up e + JfcCi/FwfJ 

Lwi^n) = y [(Co + Cpjy) Up — (Ci/C + do) \Uj’\Up e 

+{Va — l)C ( j(./J, + 4 cK a (U p — D n ) 

2 

[sign(f/ r ) ((/o + QC 0 + C MV )U^ - QKC\ |U T |U Pe 
- —cK qM (q ~V n ) 

-A^bscK Mfi (Up - z an ) - A A b±cK Mh (U p - Z bn ) 

where 

Up e = Up n - X n - Y n 

X n = X„_ie- 014 ^ As - 0.3A {/p„e-° 07 ^ As 
Yn = y n -ie" 0 - 53 ^ As - 0.7AC/p n e-- 265 ^ 2As 


(2.79) 

(2.80) 
(2.81) 


(2.75) 

(2.76) 

(2.77) 


(2.78) 


(2.73) 

(2.74) 
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with 


As ~ [A ipx - fi (cos - cos ^ n -i)] 
A 2p = ^ n ~ ^ n _ x 


A Up n = Up n - U Pn _, 


Up 



n 


D n = D n . ie - A ^ cM H K “ + AUp^-^l^nKa 

* * * 

A Up n = u Pn - Up n _ x 

V n = F n _ie-° - 5 ^ As + A g n e-°- 25 0 2As 
Aq n = qn ~ qn-l 

q n = 2 T} r (0 + <fi -f w 4- /3p + //cos ^n^tw) 


„ r = [ ~( e d + c / 2 ) 

l ~ e d 


if Uj 1 > 0 
otherwise 


V qn = V qn ^ie~ A 'l’l cM P K ‘iM + A g ne ~A^/2cM H K qM 


* 

9 = 


dq 

d® 


Aq n = q 

Z an = Z an ^ie~ A ^/ cMHK Mhh -f- AtJp n e~ A ^l 2 c M H KMhb 3 

Zf )n = Zf )n _ie A ^ ) ^ cMHKMhf}4 + A £/p e ~Aip/2cM H K Mh b4 

IC= (}_±ZM' 

!n=fn-D fn 

Dfn = D fn -ie~ As ' T f + (f n - /;_,) e-W, 

5 = &o + &i(l - fn) + sin tt /" 2 

C /V = Cat l „ - ®P n 

^ = D Pn _ 1 e~ As / T P + (C* tB - C^) e"^ 


i> 


= | (//jO-S^C^-CjvJ 


if C/V > 
otherwise 


(2.82) 

(2.83) 

(2.84) 

(2.85) 

( 2 . 86 ) 

(2.87) 

( 2 . 88 ) 

(2.89) 

(2.90) 

(2.91) 

(2.92) 

(2.93) 

(2.94) 

(2.95) 
(2-96) 

(2.97) 

(2.98) 

(2.99) 

( 2 . 100 ) 

( 2 . 101 ) 

( 2 . 102 ) 

(2.103) 
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and for C' N > C Nl 


C NV n = CNv^e-W* + (C Vn - C Vn _ i)e~ As / 2Tv 


C Vn~~ Cn l n 1 



2l 


(2.104) 

(2.105) 


with Cpy and C^v given by equations 2.71 and 2.72. The parameters JC,Q, and <t> dictate the 
nonlinear aerodynamic behavior, while the dynamic stall process is governed by C NV . For quasi- 
steady aerodynamics, the recursive parameters X n ,Y n ,U P - D„,V n ,V qn , Z an , and Z 6 „ are set to 
zero. 


Since the blade section airloads are calculated in the blade undeformed coordinate system, a trans- 
formation is required to convert these to the undeformed frame to be consistent with the virtual work 
expression, equation 2.9. For this, the aerodynamic forces in the undeformed frame are obtained from 
the coordinate transformation 

L U ) T f Lu ) 

L \ I = I - v I (2.106) 

K } I Lw J 

These force components, together with M ^ , are the motion-dependent airloads acting on a blade section 
and contribute to the virtual work expression in equation 2.9. 


2.4 Rotor Wake Modeling 


Calculation of the rotor nonuniform induced inflow is one of the most crucial steps in the prediction 
of blade aeroelastic response, blade loads, vibration, noise, and rotor performance. Currently, there are 
many inflow models available. These models try to capture, each with a different level of approximation, 
the complex structure of the vortex wake which is the most dominant feature of the inflow calculation. 
Complex models can yield accurate inflow prediction but require extensive computation. Simple models 
provide quick estimate, but somewhat less accurate results. In general, the balance between accuracy 
and computation efficiency dictates the selection of the wake model. This section describes the inflow 
models used in the present analysis. 

As a consequence of the bound circulation variation associated with the lift on each blade, vortex 
elements leave the blade trailing edge and spiral downward below the rotor disk plane. This is in 
accordance with the law of conservation of vorticity, which requires the formation of vorticity associated 
with the spanwise as well as timewise variation in the bound circulation. The spanwise variation of 
circulation produces the trailed vorticity, and the timewise variation results in the shed vorticity oriented 
radially in the wake. Such effects have been depicted schematically, as shown in figure 2.6 (adapted 
from ref. 90). The effects of shed vorticity have already been included in the unsteady aerodynamic 
model (sec. 2.3), and, therefore, this section presents methods for computing the rotor induced inflow 
associated with the trailed vortex system. 
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In the analysis, different inflow models with varying levels of complexity were included The 
simplest model assumes a uniform inflow distribution over the rotor disk and was first proposed by 
auert (ref. 91). This is an extension of the classical momentum theory to forward flight. In this 

model, the total rotor inflow is expressed in terms of the rotor thrust coefficient, advance ratio and rotor 
shaft tilt as 


A = n tan a s + X i 

C t /2 

= fJ- tan ce s + - -±1., „ 

yj \ 2 fi 2 


In another model, the rotor induced inflow is assumed to vary linearly across the rotor disk and is 
expressed as 


where 




\/A 2 + ^ 


+ k x x cos ^ -f K y x sin #) 


(2.108) 



(i - 1 V) 

. 



fcy — — 2fl 



A 


This is known as the Drees inflow model (ref. 91). 


(2.109) 

(2.110) 


Since the induced inflow is, in fact, highly nonuniform across the rotor disk, the following model 
uses a vortex ring representation. In this model, the helical path of the concentrated tip vortex for one 
complete rotor revolution is idealized as a vortex ring. Since the velocity flow field induced by a vortex 
ring can be computed analytically (ref. 92), the rotor trailed wake structures can be represented by a 
system of vortex rings cascaded below the rotor disk, as shown in figure 2.7 (adapted from refs. 90 
and 93). The resulting induced velocity at any point in the flow field is the net sum of the contribution 
from the individual vortex rings. The relative location between two rings can be determined by the 
advance ratio and the mean rotor inflow. While this model can be classified as a rigid wake model, the 
near singularities in inflow which are contributed by the tip vortices are well represented. The trailed 
vortex strength is assumed to be concentrated near the blade tip, and the strength of the vortex ring T 0 
can be expressed in terms of the rotor thrust coefficient as 


r 0 = 2flRc < X 

a 


( 2 . 111 ) 


To cover the singularity condition at the vortex core, T 0 is further modified as (ref. 91) 


r (p) = r 0 



( 2 . 112 ) 
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where p is the ratio of the radial distance from the core center to the core radius. Typically, core 
radius of 5 to 7 percent of chord length is employed. This inflow model, however, fails to capture the 
nonsymmetric characteristic of the induced inflow across the rotor longitudinal axis. 

To have a more realistic representation of the wake structures and to improve the accuracy of the 
inflow prediction, the wake analysis of CAMRAD has been incorporated into the aeroelastic analysis. 
These models capture, to various levels of approximation, the effects of the interlocking helical vortices 
of the wake and include prescribed or free wake models. The wake models consist of a set of discrete 
vortex elements representing the wake structure. The velocity induced by a vortex element is calculated 
by performing integrations over that element in the wake using the Biot-Savart law, and the resultant 
induced velocity is the sum of the contributions from all these vortex elements. The basic vortex 
elements employed by these models are the vortex line with constant or linear circulation distribution 
and the vortex sheet with uniform density or linear circulation variation. The strength of the vortex 
elements is determined by the radial and azimuthal variation of the bound circulation. 

Basically, the wake geometry consists of discrete and skewed spiral wake, one behind each blade, 
and is determined by the locations of the wake element. These wake element locations are computed 
based on the blade motions, the free-stream convection velocity, and the convection and distortion due 
to the wake self-induced velocity. In essence, classification of the wake models depends largely on the 
method used in computing the wake distortion. For the rigid wake geometry, the wake elements are 
convected downward by the mean induced inflow, and the self-distorted wake geometry is neglected. 
To include the effects of wake distortion, the prescribed wake models include wake contraction using 
assumed model or empirical data. For the free wake models, the wake induced velocity is included in 
the computation of the wake self-distorted geometry. 

Concentration of lift and circulation near the blade tip region results in the formation of the 
concentrated tip vortex which is the flow structure dominating the wake calculation. The concentrated 
tip vortex is trailed as a helix and forms behind each blade. In the wake analysis, this is modeled by 
a series of straight connecting finite-length vortices with a small core radius extended over the whole 
wake. There is also an inboard vortex sheet trailed from each blade which has a combined strength equal 
and opposite to that of the tip vortex. The inboard trailed vorticity is modeled by either a rectangular 
vortex sheet or a lattice of discrete straight-line vortex segments. 

These wake models are composed of three parts: the near wake, the rolling-up wake, and the far 
wake. The near wake is part of the wake region just behind each reference blade where the induced 
velocity is calculated. Since both the radial and azimuthal variations of the wake vorticity are important 
in the calculation of the induced inflow, a full vortex panel representation is used for the near wake. In 
the near wake region, the tip vortex is partially rolled up, and hence a line vortex with strength equal to 
a fraction of the bound circulation near the blade tip is used. The wake region just downstream of the 
near wake is the rolling-up wake. The roll-up process may not be complete as the tip vortex encounters 
the following blade. Thus, modeling of the roll-up process is required to yield the proper magnitude 
of the induced velocity, due to such an interaction. For the rolling-up wake, two wake panels with 
linear circulation variation are used. In the far wake, the roll-up process is complete, and the tip vortex 
strength is assumed to be equal to the maximum bound circulation. This wake region extends from 
the roll-up wake boundary to the far wake. For the far wake, a single vortex sheet panel with linear 
circulation variation is employed. Typically, the near wake extends from 15 to 30 degrees azimuth 
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behind a reference blade, while the rolling-up wake extends to about 30 degrees azimuth behind the 
near wake boundary. 

For a close blade vortex interaction, the induced velocity calculated using a lifting line model is 
not accurate since the induced loading would vary rapidly along the blade chord. Efficient modeling of 
this behavior is made using a correction factor based on the lifting surface solution in conjunction with 
the lifting line model. Furthermore, it has been established through experiments that the vortex core 
bursts and increases in size after such interaction. Even though the exact cause of such a phenomenon 
is still speculative, these wake models allow a vortex core to increase in size to roughly 10 percent of 
chord length after a blade-vortex interaction. 


For the free wake model, the tip vortex is modeled using straight line segments, while rectangular 
sheets or line segments are used for the inboard wake. The rigid or prescribed wake geometry is used 
for the inboard vorticity, and only the distorted geometry of the tip vortex is computed. Basically, the 
distorted wake geometry is calculated by integrating the induced velocity at each wake element. For 
a given blade azimuth angle 'F, the tip vortex geometry is given in terms of the displacement and the 
induced velocity of each wake element at the corresponding age. As the wake age is increased by AT 7 , 
the wake distortion is updated based on the induced velocities calculated at these wake elements. 

Coupling of the wake models obtained from CAMRAD and the rotor aeroelastic analysis is made 
by transferring the blade bound circulation, the blade motions, and the induced inflow between these 
two analyses. The process of finding consistent solutions for the blade unsteady airloads and the rotor 
induced inflow is solved iteratively. To start the iteration process, the blade airloads, bound circulation, 
and blade responses are computed using a uniform or linear inflow model. With the bound circulation 
and blade motions given, the wake vorticities and geometry are calculated. These are then used to 
recompute the nonuniform inflow, and the entire process is repeated until convergence is obtained. 

2.5 Blade and Rotor Hub Loads 


The resulting blade loads are computed using the force summation method. In addition to the 
distributed airloads, this method also requires information about the blade inertial loads, which can be 
determined from the blade motions. From Newton’s Second Law, the sectional inertial forces acting at 
a blade section are 

F 7 = J fp s adr]d^ (2.113) 

with F 7 = \_L^, L^, L{p L] l . L 7 , and L 7 , are, respectively, the distributed inertial forces acting in 
the x , y, and z directions in the blade undeformed frame; p s is the blade structural density; a is the 
blade acceleration relative to an inertial frame. In the blade undeformed frame, the blade acceleration 
has the form 


a — a x i + ciyj T - cl z k 

= r + Ox(flxr) + 2(nxr) (2.114) 

where r, f, and r are, respectively, the blade acceleration, velocity, and displacement vectors in the 
rotating blade frame, and 0 is the rotor rotational speed. 
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The inertial component for the blade pitching moment about the deformed elastic axis is 


M l = / fp[a>y (zi ~w)~a z ( yi v)] dr)d( 


(2.115) 


Using the sectional constants defined in section 2.2, the inertial load components can be expressed 
as (ref. 94) 


L u 171 


Ly = —m 


1 — x — u — 2v + /3 p w + e g cos 0\ (V — 

+e 5 sin (w' + (3 p ~w' + 2<?i)] 

v - v + 2u - 2 p p w - e g cos 6\ (l + 2v'^j 
—e g sin^i ($i + 2 

L J W — —m jw) + (3p (x 4- it) + 2 P p ve g cos (#i - /3p?/) 

—e g sin Q\(dp (w r + 2#i)] 

~iti — e g sin 0\ (i) — v) + e g cos 0\ (w + xfip'j 

T {j^rrii ~ ^ 77 ^ 2 ) ^OS $1 sin 0\ ^1 + 2v'^j 
+ [k^ ni sin 2 e x + /c^ 2 cos 2 O^j 2 it;'] 

The resultant distributed loads consists of the aerodynamic and inertial components; that is 


(2.116) 

(2.117) 

(2.118) 
(2.119) 


( 2 . 120 ) 


L U ~ Ly + Ly 
Ly — L A + Ly 
Lw — L^ + L\y 

M(j) — 

The resultant blade shears computed at a radial station r Q are 



dr 


and the resultant blade bending moments are 

—L v (w — w Q ) + L w {y — v 0 ) + AtiiMQ 
Lu(w w 0 ) — L w (r — r 0 + u — u 0 ) + A v' 
L u (v v 0 ) ~f Ly(v — Tq fu — u 0 ) -\- Aw' 


' M x ) 

rR 1 

My > 

= / 

1 M z J 

Jr Q j 


dr 


( 2 . 121 ) 


( 2 . 122 ) 
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where u 0 ,v 0 ,w 0 are, respectively, the blade axial, lead-lag, and flap displacements at the radial station 
r Q , and 


with 


A*n - 1 


Av /2 + Aw' 2 


2 


A /2 ' 2 / 2 
Av = v — v 0 

Ait; = w — w 0 
Av' — v — v 0 
Aw = w' — w' 0 

Similarly, v 0 and w 0 are, respectively, the blade slopes in the lead-lag and flap directions computed at 
the radial station r Q . 

At the blade root, r Q is zero, and so are the blade displacements. Therefore, the expressions for 
the blade bending moment simplify to 


r ■ 

/•is 

My 

= / n 

M z J 



-L v w + L w v + (l - M , 0 ' 

L u w — L w (r + u) + v' ' 
—L u v + L v (r + u) + w'Mfi 


dr 


(2.123) 


Furthermore, these blade loads can be expanded using the Fourier series. For instance, the radial 
shear expression for the mth blade can be expressed as 

oo 

= F x 0 + E (F™ c cosAV m + F^sinA'Pm) (2.124) 

n—1 

The Fourier series expansion of the other blade load components have a similar form. 

The rotor hub loads are obtained by summing contributions from the individual blades. In the 
rotating frame, the hub loads are expressed as 


„ ^ b 

= e G 

171 = 1 

N b 

Fy R W = E G 

m= 1 

„ N <> 

= e (*r) 

m=l 

jv 6 

= E G 

m=l 


cos A^ m - T 1 ™ sin A$ m J 

(2.125) 

sin A$ m + F™ cos A# m ) 

(2.126) 

) 

(2.127) 

1 cos AF m — My 1 sin A^ m ) 

(2.128) 
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(2.129) 


N b 


M y^)= E {M™ sin A4> m + M™ cos A\P m ) 

m = 1 

„ N » 

= e (*o 

m — 1 


where the azimuth for blade m is 
with 


= * + A$ m 
2 tt 


A$ 


m 


m- 




Substituting the Fourier series expansions of the individual blade loads into the hub load 
and using the properties of the sum of harmonics (ref. 91) yield 


N b 00 
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(2.130) 

(2.131) 

(2.132) 
expressions 
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(2.134) 

(2.135) 
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M?(®) = yf 

z p=i 



1)4> 


Mx (pN b - l)c + M V(pN b -l)s 


+ sin(pN b - 1)^ 


M 7 




(pN b -l)s ^(pA^-lJc 
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Nh 00 

+ -7T EE \ cos (P N b + 1)* 
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+ sm(pN b + 1)^ 
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(2.136) 
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(2.137) 

(2.138) 


For the rotating inplane hub shears, it is worth noting that 


and 


F R 

= F R 

x (pN b -l)c 

y(pN b ~i) S 

F R 

= -f r 

x {pN b -l) s 

y(pN b ~ i)c 

F R 

= —F r 

x (pN b + \)c 

y( P N b +i)s 

F R 

= F, r 

x (pN b + l)s 

y(pN b +i) C 


Therefore, the rotating inplane hub shears are uniquely determined by either F R or F R alone. 
Similar properties exist for the rotating inplane hub moments M R and M R . Also, the rotating hub 
loads occur only at certain frequencies. The frequency contents of the vertical hub shear F r and those 
of the torque M R are the harmonics at pN b per rev, where p is an integer, and p = 0, 1, 2, 3, . . ., whereas 
the rotating inplane components occur at (pN b — 1) and (pN b + 1) per rev, where p is a strictly positive 
integer. Therefore, only certain harmonics of the blade loads are summed up while others cancel at the 


44 



rotor hub. The filtering of the blade load harmonics at the hub is caused by the equal phasing differences 
of A$ m between individual rotor blades. Furthermore, it has been assumed that the rotor is tracked. 
For a tracked rotor, all the blades exhibit exactly similar behavior except for the phase differences 
of A4> m . 

Similarly, the expressions for the rotor hub loads in the fixed frame are 


F X 0*0 = E i^x 1 cos — FV 1 sin ® m ) (2.139) 

m = 1 
N b 

F Y 0*0 = 12 ( F ? sin(& m + i^cos'i' m ) (2.140) 

m — 1 

„ N » 

F Z^) = E (FT) (2.141) 

m — 1 
„ N h 

= (M^cos^ m -M^sin^ m ) (2.142) 

m — 1 
__ X b 

m y(^)= ]£ (M™ sin \I/ m + M™ cos (2.143) 

ra= 1 
„ N b 

**#(*) =E(*C) (2.144) 

m = 1 


For a tracked rotor, the Fourier series expansions for these hub loads can be expressed as 

A T 00 r ^ 

W = T ^ l C0SpAr ^ F HvN b - l)c + 

P=1 L J 
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FyW = Y E {cospTV^ [-F x ( + F„ ( 

+ sin P Ar 6 $[F X(pW6 _ 1)c + F y(pJV6 _ 1)s ]} 

AT, 00 r r -1 

+ ~2~ F x {pN + i )s + F y( pN + i) c 
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+ sin PN b * [^ (pW6 _ 1)c +M y(pWfc _ 1)s ]} 

+ T E : {cospiVi* ^ s + M„ (p 
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+ sin pN b ^i [~M X(pNb+l)c + M y(pW6+1)s ]} (2.148) 

whereas the expressions for the vertical shear and the rotor shaft torque in the fixed frame are 
the same as those in the rotating frame. 

As these hub loads are transmitted to the fuselage, the steady components are the rotor thrust, 
longitudinal and side forces, rolling and pitching moments, and the rotor shaft torque. The higher 
harmonics contribute to helicopter vibration. In particular, helicopter vibration occurs at pN ^ per rev, 
where p is a positive integer. As in the case of the rotating hub loads, only certain harmonics of the 
blade loads are transmitted to the hub fixed frame. Again, this happens because of the equal phasing 
between individual rotor blades. Furthermore, since the load harmonics tend to lower magnitudes at 
higher frequencies, helicopter vibration occurs primarily at per rev. 
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2.6 Coupled Trim Analysis 


Typically, the vehicle trim controls and the rotor steady responses are solved independently. The 
solution method is usually performed in the following manner: (1) For a prescribed flight condition, 
the trim controls are computed based on a rigid flapping blade assumption. For this, the rotor induced 
velocity is calculated based on the momentum value, and the quasi-steady approximation is used to 
calculate for the blade loads. This would lead to a crude approximation to the actual trim solutions. 
(2) Based on these trim solutions, the blade steady responses are then calculated. For this step, the 
elastic blade model together with more refined wake and aerodynamic analyses are employed. Then the 
rotor hub loads and a new vehicle equilibrium condition are computed. The new equilibrium condition 
is usually quite different from the one prescribed initially. 

For a given flight condition, the trim controls and the blade steady responses must satisfy both the 
blade and the vehicle equilibrium conditions. The method incorporating the process of simultaneously 
solving for the blade responses and the trim controls is referred to as the coupled trim analysis. In this 
analysis, the blade steady responses are solved iteratively. With each new blade response update, the 
vehicle would attain a new equilibrium condition. Therefore, the trim controls are also updated so' that 
the overall solutions satisfy the prescribed vehicle equilibrium conditions. 

In this analysis, an uncoupled trim solution based on the rigid flapping blade assumption (ref. 95) 
is used as an initial guess for the coupled trim analysis. With these trim controls, the blade steady 
responses are calculated. Using the blade responses, the rotor hub loads and a new vehicle equilibrium 
condition are computed. The trim controls are then updated based on the new vehicle equilibrium 
condition. The process is repeated until the overall solutions converge. For a given flight condition, 
the converged trim control and blade response solution satisfy simultaneously the overall vehicle steady 
force and moment equations. 

In particular, the general expression for the vehicle equilibrium condition can be expressed as 

F = 0 (2.149) 

where the exact form of — [F\ , . . . , F n J depends on the types of trim conditions considered. In the 
analysis, two types of trim are considered: propulsive and wind tunnel. The propulsive trim simulates a 
helicopter in a steady-level flight condition as shown in figure 2.3, and the vehicle equilibrium conditions 


must satisfy the following force and moment equations 

F l = F X 0 + D F cosa HP - T sin a s (2.150) 

f 2 = Fj? + Y F cos (p s +T sin <f) s (2.151) 

*3 = F z 0 ~ T cosa s cos(f) s - D F cos <fi s sin a jjp + Y F sin (f) s (2.152) 

f 4 = + Mx f + T (h cos a s sin (fi s - y CG cos 

+ Y f (h cos a s cos <j) s + y C G sin <t> s ) (2. 1 53) 
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^5 = My 0 + M Vf +T (h sin a s — xqq cos a s ^j 

+ D F [~h cos a HP + x CG sin a HP ^j (2.1 54) 

where F\, F 2 , and F 3 are, respectively, the vehicle force equilibrium equations in the X , Y and Z 
directions in the hub-fixed frame, and F 4 and F 5 are the vehicle rolling and pitching moment equilibrium 
equations, respectively. The yawing moment equilibrium equation is not considered since the tail rotor 
dynamics are neglected in this analysis. Also, D F is the fuselage drag; Y F is the fuselage side force; X 
is the rotor thrust; x CG , y CG , and h are, respectively, the relative location of the rotor hub center with 
respect to the vehicle center of gravity in the X , Y, and Z directions; a s and <j> s are the longitudinal 
and lateral shaft tilts, respectively. Furthermore 


a s — <*HP ~ &FP 

where a HP is the hub plane tilt relative to the flight direction, and 0 FP is the flightpath angle relative to 
an axis perpendicular to the weight vector. For propulsive trim, the unknown quantities to be determined 
from the vehicle equilibrium equations are 

u = [a Si 0s, 0.75, 0i c ,6\ s \ (2.155) 

where the rotor trim parameters 0 75 , 0 lc , and 6>i s govern the blade pitch schedule. With no higher 
harmonic control, the blade pitch schedule can be expressed as 

0(x,^>) — 6(75 +$t w (x -0.75) -f- 0\q cos + 6i s sin ^ (2.156) 

where x is the nondimensional radial length; 0t w is the blade linear pretwist; # 75 is the collective pitch 
angle at 75 percent radius; 0\ c and 9\ s are, respectively, the lateral and longitudinal cyclic pitch angles. 

Wind tunnel trim simulates the rotor testing conditions in a wind tunnel. For such conditions, the 
longitudinal and lateral shaft tilt angles are typically prescribed. Therefore, thp longitudinal and lateral 
force equilibrium equations need not be satisfied. As such, the unknown vector to be determined from 
the remaining equations F 3 , F 4 , and F5 is 

« T = (0.75. 01c, #ls J 


For certain wind tunnel tests, the trim procedure sometimes requires the suppression of the cyclic 

flapping angles measured at a reference blade station. For such cases, the expressions for Fa and Fk 
become 


d 1 r p 

F 4 = wF. = - / w cos ^ d^J 

7 r Jo 

j-, d 1 r p 

F5 = u>i s — — w n sin \1> 

7 r Jo 


(2.157) 

(2.158) 


Finally, the solution technique is based on a modified Newton method, and, as such, the trim 
controls are updated at the ith iteration as 


u i+ 1 = u i + A11* 


(2.159) 
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where 


F K) (2.160) 

u=u 0 

where u 0 are the trim controls obtained initially using the rigid flapping blade solution, and the Jacobian 
dF/dn is calculated using the finite-difference approach. Due to computational requirement, the Jacobian 
for trim is computed only once initially and is not updated for subsequent iterations. 

To improve numerical convergence, a numerical damping factor R is introduced into equation 2.159 
which can be expressed as 

u i+l = + R All; (2.161) 

Currently, there are two forms for R. In the first form, R is given by 

R — 1 ~~ r 

where 0 < r < 1. Additionally, R can also take the form 



Typical value of ( r is l/y/2. 


R= 1 - exp(-Cr^) 


2.7 Finite Element Discretization in Space 


This section describes the finite element method used in the derivation of the blade governing 
equations. In particular, the finite element analysis is used here to remove the spatial dependence from 
the blade governing equations. For this, the rotor blade is discretized into a number of beam elements, 
and shape functions are used to relate the elemental degrees of freedom. 

From the Hamiltonian formulation, the integrand of equation 2.3 can be expressed as 


A = 8U-8T~8W 


(2.162) 


and in discretized form 

n 

A = E A i (2-163) 

1=1 

where 


A i = 8Ui - 8 % - 8Wi (2. 1 64) 

and 8Uj, 871 , and are, respectively, the variations in strain and kinetic energies and virtual work 
contribution of the ith element, and n is the number of beam elements. As shown in figure 2.8, each 
beam element has 15 degrees of freedom and consists of two end nodes and three internal nodes. Each 
end node has 6 degrees of freedom which are: u, v, v f , w, w f , and </>. There are two internal nodes for 
the u degree of freedom and one internal node for the <j> degree of freedom. There are two main reasons 
for choosing 0 over (f) as the torsional degree of freedom: (1) the choice of </> preserves the banded 
structure of the global finite element matrices since the integral term v n w r dx is implicitly included 
in the formulation, and (2) the aerodynamic loads depend explicitly on the geometric twist </>. At the 
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element boundary, there is a continuity of displacement for the axial, lag, flap deflections, and twist 
and a continuity of slope for the flap and lag bending. 

In the formulation, the deflections at any point within an element can be expressed in terms of the 
shape functions and the nodal degrees of freedom. For the ith element 


u(s ) 
v(s) 
w(s) 

, k s ) 


> = H(s) q e 

> 


where the nodal vector q e consists of the 15 degrees of freedom for the beam element, and 


T 

q e = \u\ , u 2 , u 3 ,U 4 ,Vi , v [ , v 2 , v' 2l wi , w [ , iu 2 , 4> 3 J 


(2.165) 


(2.166) 


The shape function matrix H consists of the shape functions for the beam element degrees of 
freedom and has the form 


H = 


H%(s) 0 » 

0 Hi (.,) 0 

0 0 H'l(s) 

0 0 0 


0 

0 

0 

H l (*) 


where the shape function vectors are 


H u (s) - { 


\ —4.5s 3 + 95 2 — 5.5s + 1 ) 
13.5 s s - 22.5s 2 + 9s 
-13.5s 3 + 18s 2 -4.5s 
4.5s 3 — 4.5s 2 + s 


H v (s) = H w (s) = 


{ 2s 3 - 3s 2 + 1 
^i(s 3 2 s 2 -j~ s) 

-2s 3 + 3s 2 
k(s 3 - s 2 ) 





(2.167) 


(2.168) 


(2.169) 


(2.170) 


The local spatial coordinate s of the ith element is given by 


Xi 


where is the length of the ith beam element; x % is the local coordinate of the beam element, and 
0 < The shape functions for the lag and the flap bending deflections are the Hermitian 
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polynomials which allow continuity of deflections and slopes across the beam element. For the axial 
deflections and elastic twist, Lagrangian polynomials are used as shape functions which allow continuity 
of displacements across the beam element. 


Similarly, the virtual displacements 8u , 8v, Sw, and 8cj> over the zth element can be expressed in 
terms of the same set of shape functions as 


j M 5 ) 

I 8v(s) 
8w(s) 

k 6 4 >( s ) 


= H(s) <5q e 


(2.171) 


where the virtual nodal vector 8q e is the variation of q e . Substitution of equations 2.165 and 2.171 into 
equation 2.164 gives 

A { = Ai(q e ,q e ,qe,<5qe) (2.172) 

Following equation 2.163, the assembly of the beam elements yields a set of nonlinear ordinary 
differential equations in terms of the nodal displacements. Let q denote the global displacements for 
the rotor blade and <5q the associated global virtual displacement; the global equations can be expressed 
implicitly as 

A = A(q, q, q, <5q) (2.173) 

Derivation of the blade governing equations is completed by applying the geometric boundary 
conditions to equation 2.173. The rotor configuration, such as hingeless or articulated, dictates the 
displacement boundary conditions near the blade root. In particular, u, v\ w , w f , and </> are zero at the 
blade root for a hingeless rotor. For an articulated rotor with the flap and lag hinges coincident, u, v , w, 
and <fi are zero at the hinge. In the analysis, further modifications have been made to include the pitch 
link stiffness, different offsets for the flap and lag hinges, and the pitch link location. 

After applying the boundary conditions, equation 2.173 takes the form 

A = 6q T (M G q + C G q + K G q - F G ) (2.174) 

where Mq,Cg,Kg are the global mass, damping, and stiffness matrices, respectively, and Fg is the 
global force vector. Substituting the above equation into equation 2.3 (note that the blade steady response 
is periodic for a rotor operating in steady-level flights), these give 

f2lT 

J Q <$q(M G q + C G q + K G q - F G ) d$! = 0 (2.175) 

Note that the above equation has been expressed in terms of the rotor temporal coordinate Tf. 

2.8 Solution Procedure 


This section describes the solution methods for the blade response and the vehicle trim controls. In 
the analysis, the blade governing equation (eq. 2.164) is transformed to the normal mode space using a 
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few coupled natural vibration modes. This is performed to reduce the computational time. The resulting 
normal mode equations are then solved using a finite element in time method. With the temporal finite 
element method, the blade normal mode equations are transformed into a set of nonlinear algebraic 
equations which are solved using a Newton method. 

Starting with the description of the nondimensional procedure, this section presents the normal 
mode transformation and the temporal finite element method used in the determination of the blade 

aeroelastic responses. Finally, implementation of the aerodynamic model into the aeroelastic analysis is 
also discussed. 

2.8.1 Nondimensional Procedure 


The parameters used in the analysis are normalized according to the following procedure 


Quantities Nondimensional parameters 


Length 

R 

Time 

I/O 

Mass/length 

m 0 

Velocity 

QR 

Acceleration 

Q 2 R 

Force 

rrioVlP'R? 

Moment 

m 0 f2 2 i? 3 

Energy or work 

m 0 n 2 R 4 


where m 0 is a reference blade mass per unit length. 

2.8.2 Normal Mode Equations 


To reduce computational time, the blade finite element equations are transformed into the normal 
mode space using a few coupled natural vibration modes. The blade natural vibration modes are obtained 
by solving an algebraic eigenvalue problem. By neglecting the external loads and the damping matrix, 
the blade governing equation can be expressed as 

M G<1 + K G<I = <> (2.176) 

Let 


Then equation 2.176 can be arranged to yield the following algebraic eigenvalue problem 

w 2 q = (Mfj ) “ 1 K?;q (2. 1 77) 

where the superscript s denotes structural quantities. These matrices are the structural inertial and 
stiffness matrices for the rotor blade and are symmetric. This eigenvalue problem can be solved using 
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the Jacobi method. The resulting eigenvalues, being are real, and the associated eigenvectors are 
also real and orthogonal. Physically, the eigenvectors represent the blade natural modes, and the square 
root of the associated eigenvalues represents the blade natural frequencies. 

To apply the normal mode transformation, the blade global displacement vector q is related to the 
modal displacement vector p by the normal mode transformation 

q = $P (2.178) 

where $ is an Nq x m matrix formed with the m blade free vibration modes. Applying the normal 
mode transformation to equation 2.175 results in 

r2ir 

<5p (Mp -f Cp + Kp — F) d# = 0 (2.179) 

V U 

where 


M = M(j3> 

c = $ T C G $ 

K = •J> 2 K G 'I> 

F = $ t F g (2.180) 

are, respectively, the modal mass, damping, stiffness matrices, and load vector. In general, M, C, 
and K contain periodically time- varying coefficients; that is 

M(tf) = M(\P + 2n) 

C(W) - C(tt 4- 2tt) 

K(tf ) = K (tt + 2 t r) 

and F is nonlinear. For periodic solutions, F should also be periodic. 

2.8.3 Finite Element Discretization in Time 


The blade modal equation (eq. 2.179) can be integrated by parts to yield 


( ip V T < F - Cp - Kp 
Jo [ 6p j | Mp 




(2.181) 


For a periodic solution, the right hand side of the above equation is zero, and the above equation 
can be expressed as 



(2.182) 
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where 


and 



O — / F - Cp - Kp 1 
1 Mp / 


Using the finite element method, the time interval for one rotor revolution of 2n is divided into a 
number of time elements in a manner similar to the spatial finite element method described in section 2.7. 
This is illustrated in figure 2.9. The blade governing equations can then be expressed as the sum of the 
elemental expressions for all the time elements. In particular, equation 2. 1 82 takes the form 



(2.183) 


where - 0, tfjVe+l = 2 tt, and N e is the number of time elements being used. The first order Taylor 
series expansion of equation 2.183 about a steady state value of y Q = |pj pTjT i s 



SyfQi(y 0 + Ay) d V 


4-A fVi+l T 

z2 L. 6 y i [Qi(yo) + K £i (y G )Ay] dV = 0 

2—1 


(2.184) 


where 


Kti = 


r <9F 
<5p 


-K 

0 



(2.185) 


For the zth time element, the temporal variation of the modal displacement vector can be expressed 
in terms of the shape functions and the temporal nodal displacement vector & as 

PiW = (2.186) 

where the local temporal coordinate s for the zth time element is 

# - Va 

s = — — 

*i + 1 - *i 

and 0 < 5 < 1. Note that ^ i+i - is the time span of the ith time element. In equation 2.186, H(s) 
is the temporal shape function matrix which has the form 

H . . . , H nt + 1 1777 , J 

where I m is the m x m identity matrix, and m is the dimension of the modal displacement vector. In 
the above expression, the subscript n £ refers to the order of the polynomial used in approximating the 
temporal variation of the modal displacement vector. In particular, if an nth order polynomial is used 
in the approximation, then n + 1 nodes per degree of freedom are required to completely describe the 
variation of p^W) within the time element. Therefore, H is an m x m(n £ + 1) matrix, and is an 
m(m 4-1) vector. 
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In this analysis, fifth-order polynomials are used in the approximation, and this results in the 
requirement of six nodes per degree of freedom for each time element. Between elements there is con- 
tinuity of displacements. With this type of approximation, the shape functions used are the Lagrangian 
polynomial and are expressed as 

Hi(s) = (-625 s 5 + 1875s 4 - 2125s 3 + 1125s 2 - 274s + 24) /24 
H 2 (s) = (3125s 5 - 8750s 4 + 8875s 3 - 3850s 2 + 600s) /24 
H 3 (s) = (—3125s 5 + 8125s 4 - 7375s 3 + 2675s 2 - 300s) /12 
H 4 (s) = (3125s 5 - 7500s 4 + 6125s 3 - 1950s 2 + 200s) /12 
#5(s) = (-3125s 5 + 6875s 4 - 5125s 3 + 1525s 2 - 150s) /24 
H 6 (s) = (625s 5 - 1250s 4 + 875s 3 - 250s 2 + 24s) /24 

With this description, the modal velocity vector can be expressed as 

Pi = H($) & (2.187) 


In a similar fashion, the expressions for the virtual modal displacement and velocity vector are 

£Pi = !*(*)$£,• (2.188) 

<5pi = H(\P)£& (2.189) 


Substituting equations 2.186-2.189 into equation 2.184 yields 


where 



(2.190) 


(2.191) 


For periodic solution, the boundary condition for the temporal finite element equation is 


p(0) = p(2rr) 


Finally, since the <5£; are arbitrary for i = 1, . . . ,N e , equation 2.190 takes the form 


where 


Q G + Kp £ G = 0 




N t Q, 


(2.192) 
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and 


K f = E /^’ +1 N T K«N d<£ 

% — 1 ^ 

N e 

A^ G = Y, A ft 

2=1 

Equation 2.192 is a set of nonlinear algebraic equations and can be solved using Newton’s method. 
In particular, are solved iteratively, and the overall responses are updated using 

£ G i+i = + a 


During the iteration process, the vehicle trim controls are also updated according the procedure 
outlined m section 2.6. The converged trim control and blade steady response solutions satisfy si- 
multaneously the vehicle steady force and moment equations (i.e., eq. 2 . 6 ) and the blade equilibrium 
conditions given by equation 2.182. 

Implementation of the unsteady aerodynamic model described in section 2.3 required calculation 
of the aerodynamic states. Since the constructions of Q G and Kp of equation 2.192 are essentially 
a time-marching algorithm, the recursive formulas in the expressions for the unsteady airloads can be 
adapted to the current analysis. Also, the unsteady aerodynamic states such as X n ,Y n , or D n are 
functions of the down wash velocity Up or the airfoil pitch rate q which, in turn, are dependent on the 
blade motions and the trim controls. Therefore, as the blade responses and trim controls are updated, 
these aerodynamic states are also updated. In the limit, as the solutions converge, these aerodynamic 
states also converge to their corresponding values. 

2.9 Higher Harmonic Control Model 


In the HHC analysis, a linear frequency-domain quasi-steady model in terms of a transfer matrix 
is adopted to relate harmonics of the HHC inputs to harmonics of the vibratory hub loads. The optimal 
HHC inputs are calculated by minimizing the harmonics of the vibratory hub loads expressed in terms 
of a quadratic performance index. Using the multiblade coordinate transformation, HHC inputs in the 
rotating frame are transformed to the fixed frame in terms of swashplate motions. The higher harmonic 
swashplate motions are generated using a number of auxiliary servo-actuators, as shown in figure 2 . 10 . 
The hydraulic power for the actuators is then calculated for the specified HHC inputs. 

The higher harmonic blade controls are implemented by exciting the blade pitch at higher harmonics 

of rotational speed. For vibration control, the blade pitch schedule, expressed in the rotating frame, is 
given by 


0(x, #) = 6 ' 75 + 6 tw (x - 0.75) + #i c cos ^ + 6 > ls sin 
N b +l 

+ (&nc cos n# + 9 ns sin n^) (2.193) 

n—N b -l 
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where 0 nc and 0 ns are, respectively, the cosine and sine components of the HHC inputs. Note that 
the summation of the higher harmonic control inputs covers three frequencies: — 1,A^, and 

Nfy + 1 per rev. 

For the HHC model, a transfer matrix T is used to relate the input pitch harmonics and the output 
hub load harmonics; that is 

z = z o + T0 (2.194) 

where the response vector z has ten components and consists of sine and cosine components of the 
vibratory hub loads — the vertical, longitudinal, and lateral hub shears and the pitching and rolling 
moments, and z 0 is the uncontrolled response vector. In particular, the components of z or 2 0 are the 
higher harmonic components of the hub loads, as presented in section 2.5. The HHC input vector 0 
has six components and consists of sine and cosine components of the higher harmonic blade pitch at 
the three distinct frequencies of — 1 ,iV&, and N ^ + 1 per rev. In this analysis, the transfer matrix is 
generated using a finite-difference approach. 

Note that the HHC model as given by equation 2.194 is equivalent to the first order Taylor series 
expansion of a more general static nonlinear description which is given by 

z-f(zo,0) (2.195) 

From this, the transfer matrix T is simply the Jacobian of /, calculated about a referenced input 0 

(2.196) 


2.10 Higher Harmonic Control Algorithms 



For minimum vibration, the optimal control input 0 is determined from the minimization problem 
of the quadratic performance index J 


min J — min (^z T W z z + ^0 T Wg$\ 


(2.197) 


subject to equation 2.194. 


In the above expression, W z and W e are weighting matrices for z and 0 , respectively. The optimal 

control 0 depends on the type of controller used. Two types of controllers are used in this analysis: 
open-loop and closed-loop deterministic controllers. 


For the open-loop controller, only the global model is applicable, and the optimal control 0 depends 
only on the uncontrolled response; that is 

$ 

9= Cz 0 (2.198) 

where 


C = ~(T t W z T + W e )~ l T T W z 


(2.199) 
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and the transfer matrix is calculated about the uncontrolled response value. If the HHC system is linear, 
then the optimal control 6 would suppress completely the vibratory response in one cycle. 

For the closed-loop controller, the general relation as given by equation 2.194 is valid for each 
HHC cycle and is expressed as 


z n = Zn - 1 + T n (6 n - 9 n -\) (2.200) 

In this formulation, the dynamics of the HHC system have been neglected. This is justified by assuming 
that the rotor system is stable and that the HHC inputs are updated only after all the transient responses 
have died out. Furthermore, this assumption is justified based on the wind tunnel results of reference 33 
which showed that the time constant for the HHC system is rather small, being less than two rotor 
revolutions. With this scheme, the HHC inputs are computed based on the current response vector 

9 n = 0 n _! + (1 - r)C n z n _i (2.201) 

where 

Cn = ~ (Tn W Z T? + W e ) _1 Tj W z (2.202) 

and r is the control rate limiting factor with values ranging from 0.0 to less than 1.0. The factor (1 - r) 
constrains the HHC control rate and helps to reduce large control excursions. This is useful from the 
actuator hardware point of view, and a suitable value of r can be used to stabilize the feedback system. 

Two types of models are used to implement the closed-loop HHC system: local and global. In the 
local model, the relation in equation 2.200 is applicable for a sufficiently small range of control values, 
@n ~ &n— 1- For the global model, the HHC model is assumed to be linear for the complete range of 
control value, and thus the transfer matrix remains constant for all cycles. For both models, the transfer 
matrix is calculated initially using the finite-difference approach. For subsequent HHC cycles, the global 
model uses the same transfer matrix calculated initially, while the local model updates it using 


rp rp , ( Vn T n S n )Sn 

-‘-n+l 1 n i T 

^Ti 

(2.203) 

where 


Vn — z n+ 1 — z n 

(2.204) 

and 


s n — CnZn 

(2.205) 

This updating scheme, generally known as the secant method, was proposed by Broyden (ref. 96) 
as a solution technique for systems of nonlinear equations. Even though the transfer matrix is computed 


using the local control response, this updating method is quite different from the widely adopted Kalman 
filter used in many HHC applications. 
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2.11 Actuator Power Required for HHC Systems 


To calculate the actuator power required for an HHC system, the optimal HHC inputs are trans- 
formed in terms of swashplate motions (ref. 97) 


' r* ' 

long,c 
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-1 

0 ' 


' e (N b -l)c ' 
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1 _ 


. e (N b +l)s , 


(2.206) 


where C col , C iong , and C iat are the N b P swashplate angles in the collective, longitudinal, and lateral 
modes, respectively, and the subscript c refers to cosine component and 5 refers to the sine component. 


The actuator power P R generated by a hydraulic pump the is product of the pump pressure p e and 
the average pump flow rate Q av (ref. 98) 


Pr — Pe Qav 


(2.207) 


For HHC applications, the pump flow rate must be sufficient to excite the N b P swashplate motions 
given by the prescribed HHC inputs. For a particular actuator, the flow rate Q depends on the swashplate 
oscillatory frequency and the actuator displacement; that is 


Q = N b QS n A e 


(2.208) 


where A e is the effective piston area, and S n is the actuator displacement. The average pump flow rate 
Qav is calculated assuming that an accumulator is used to smooth out the flow; thus 


Qav — N f) S<n Af 

7 r 


(2.209) 


Throughout the helicopter flight regime, the actuators must operate at a constant fluid pressure 
selected to accommodate the maximum possible pitch-link load. Such loads occur at severe flight 
conditions such as a high g maneuvers or a high speed condition where a major part of the rotor disc 
is exposed to dynamic stall. For this analysis, it is assumed that the fluid pump pressure is a known 
constant. Therefore, the actuator power required for an TV^-bladed rotor can be expressed as 

2 

PR = - ^ Pe c n S n (2.2 10) 

7 r 

where c n is a constant dependent on the radial locations of the actuators and the pitch link, and S n 
is the accumulation index for the swashplate higher harmonic motion amplitudes. If the swashplate is 
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oscillated using four equally spaced hydraulic actuators, then S n can be expressed as 

= \A^co/,c + Cldt,c) 2 + ( C co [ s + Ci at , s ) 2 
+ \fJCcol,c - Qat,c ) 2 + ( C coUs - Ci attS ) 2 
+ yf{Ccal,c + C long,c ) 2 + (C co ^ s + Ci ong ^ s ) 2 

+v (C col , c - C long c ) 2 + (C cd , s - Ci ong s ) 2 (2.21 1) 


The actuator power for an HHC system is therefore dependent on the number of blades, the rotor 
rotational speed, the hydraulic pump pressure and diameter, the locations of the pitch links and actuators, 
and the HHC inputs. Since most of these parameters are fixed for a given rotor design, the associated 
actuator power is conveniently represented in term of a power index defined as 

p I = E r 

( 2 / 7 t) N^QAePeCn 

This gives 

Pj — S n 

Thus, for a given actuator design and pitch-link configuration, the actuator power depends only 
on the swashplate motions which are related to the HHC inputs. Also, it should be noted that this 

assumption is valid only when the maximum pitch-link loads are not significantly affected by the HHC 
inputs. 


( 2 . 212 ) 

( 2 . 213 ) 
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3 RESULTS AND DISCUSSION 


In section 3.1, the predicted vibratory hub loads are correlated with experimental data to validate the 
analysis and to set general guidelines for further refinements to the analytical model. HHC effectiveness 
in suppressing the vibratory hub loads are then evaluated. Secondary effects due to HHC are also 
investigated including HHC penalties on blade stresses, control loads, rotor performance, and influences 
of HHC on blade stall. For the cases where the simple fixed-gain controller is not effective in suppressing 
the vibratory hub loads, results using an adaptive-gain controller is presented. 

In section 3.2, a parametric study is carried out to investigate the influence of blade and rotor 
characteristics on the vibratory hub loads and the HHC actuator power. Comparative study is made 
for three different three-bladed rotor systems including an articulated rotor, a soft-inplane hingeless 
rotor, and a stiff-inplane hingeless rotor. Then, detail parametric investigations are carried out for a 
generic four-bladed soft-inplane hingeless rotor operated in a transition flight regime and a three-bladed 
articulated rotor operated in a quasi-steady maneuvering flight condition. These studies are carried out 
by varying the blade design parameters such as the blade flap, lead-lag, and torsion stiffness distributions, 
the chordwise offset of the blade center of mass from the elastic axis, the chordwise offset of the elastic 
axis from the quarter chord, and the blade linear pretwist. 

3.1 Correlation Study and HHC Effectiveness 


3.1.1 Correlation Study 


To validate the analysis, a correlation study of predicted vibratory hub loads with experimental 
data is carried out. The data were obtained from a one-sixth dynamically scaled three-bladed articulated 
rotor model of the CH-47D helicopter which was tested by Boeing Helicopter Company in their wind 
tunnel (ref. 46). In the experiment, the vibratory hub loads were measured in the rotating frame using 
a strain-gauge balance located under the rotor hub. For this rotor, the dominant vibratory hub loads 
causing airframe vibration are the 3 P vertical and the 2 P and 4 P rotating inplane hub shears. The 
vibratory hub moment components for this rotor are small, and these are not considered in the correlation 
study. 

In the validation study, six beam elements are used to model the blade, and the corresponding 
blade parameters are presented in table 3.1. The calculated rotating blade natural frequencies for the 
coupled modes used in the correlation study are presented in table 3.2. In order to apply the nonlinear 
unsteady aerodynamic model (refs. 81 and 82), static characteristics for the VR-12 airfoil are used, 
while the three time constants used to model the unsteady aerodynamic characteristics are taken from 
reference 81. The basic trim controls are calculated corresponding to the wind tunnel conditions with 
a specified advance ratio, rotor thrust, and shaft tilt angle. For wind tunnel trim calculations, the cyclic 
flapping angles measured at the flapping hinges are suppressed to zero. 



Correlation results are presented in terms of the magnitudes of the vibratory hub loads at the 
specified harmonics. For instance, the magnitude of the 3 P vertical hub shear is given by 



Figures 3.1 3.3 show the correlation of vibratory hub shears over the range of forward speeds 
from hover to 188 knots. The thrust coefficient-solidity ratio ( C T /a ) is held at 0.08. Calculated 
results are shown for three aerodynamic models: (1) quasi-steady aerodynamics and linear inflow model 
(Drees), (2) nonlinear unsteady aerodynamics (including dynamic stall) and Drees’ inflow model, and 
(3) nonlinear unsteady aerodynamics and free wake model. 

Figure 3.1 shows the correlation of the 3 P vertical hub shear. Predicted results obtained using 
nonlinear aerodynamics and free wake modeling agree quite well with measured data at low speed as 
well as high speed flight regimes. At forward speeds above 160 knots, the sharp rise in the 3 P vertical 
hub shear is attributed to dynamic stall, and this effect is reflected in the analytical results using nonlinear 
aerodynamics and Drees inflow model. The sharp peak in the 3 P vertical hub shear in the transition 
flight regime (around 40 knots) is satisfactorily predicted using the free wake analysis. It is interesting 
to note that the simple model with quasi-steady aerodynamics with Drees’ inflow correlates quite well 
with experimental data at the forward speeds above 100 knots. However, the calculated results using 
nonlinear aerodynamic and free wake models yield the best overall correlation with experimental data. 

Figure 3.2 shows the correlation of the 2 P rotating inplane hub shear. For the range of forward 
speed below 100 knots, the calculated results using all three aerodynamic models are almost identical 
except for the small peak in the 2 P inplane hub shear in the transition flight regime, captured only with 
the free wake analysis. At the forward speeds above 135 knots, the results using nonlinear aerodynamics 
with both inflow models correlate better with the experimental data than the results with quasi-steady 
aerodynamics and Drees’ inflow. 

Correlation of the 4 P rotating inplane hub shear is shown in figure 3.3. This figure also shows that 
the results obtained using the three aerodynamic models are almost identical for the forward speeds below 
100 knots. At the forward speeds above 160 knots, the results obtained using nonlinear aerodynamics 
and Drees inflow overpredict the 2 P rotating inplane hub shear, while the results with quasi-steady and 
Drees’ inflow underpredict the experimental values. At the high forward speed range, the calculated 
results using nonlinear aerodynamics and the free wake analysis yield the best correlation of the 4 P 
rotating inplane hub shear with the experimental data. Based on the results of the validation study, 
subsequent investigations were all carried out using nonlinear aerodynamics and the free wake model. 

For the flight conditions considered above, optimal closed-loop HHC is introduced in the analysis 
to suppress the vibratory hub loads. The fixed-gain closed-loop controller is used with no control rate 
limiting factor set to zero. This controller suppresses completely and simultaneously all three vibratory 
hub shear components over the forward speed range from hover to 1 60 knots. These results are shown 
in figures 3.4— 3.6 in which the predicted vibratory hub shears calculated using nonlinear aerodynamics 
and the free wake model are compared with the experimental values. Also, the results are presented 
for the cases with both HHC off and HHC on. In the experimental data, there are some residues in the 
levels of the vibratory hub load being suppressed. These residues may be attributed to the measurement 
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and process noise which are not included in the calculated results. In fact, if measurement noise is 
included in the simulation, then the HHC controller suppresses the vibratory hub loads to the noise 
level. For the case at 188 knots, the HHC amplitudes required to suppress completely the vibratory hub 
shears were very large (on the order of 5 degrees), and a converged numerical solution could not be 
achieved. It is suspected that such a large higher harmonic blade pitch might have induced a massive 
blade stall which caused the numerical solution to diverge. Such effects will be investigated in more 
detail in section 3.1.3. 

The higher harmonic swashplate motions required to suppress the vibratory hub shears are shown 
next. The HHC inputs are transformed in terms of the harmonics of the swashplate motions in the 
collective, longitudinal cyclic, and lateral cyclic modes. The magnitudes of these swashplate motions 
are then calculated, the results being presented in figure 3.7. These results indicate that, except in the 
transition flight regime, the amplitudes required by both cyclic modes are generally higher than the 
collective amplitudes. Above 100 knots, it is worth noting that both of the cyclic swashplate ampli- 
tudes increase significantly with the forward speed, while the collective component remains essentially 
constant. 

The actuator power required by the HHC system, expressed in terms of the actuator power index, 
is presented in figure 3.8. The actuator power index represents an accumulative index for the swashplate 
motions. Also, it represents approximately the combined HHC amplitudes required by all the blades 
of a rotor. The results show that the power index increases moderately from hover to a peak value in 
the transition flight regime (at 30 knots), then decreases slightly up to the forward speed of 60 knots, 
and increases again quite significantly for the high speed range (above 60 knots). The power index at 
160 knots is about three and a half times that in the transition flight regime, at 30 knots. 

Further investigations of HHC swashplate motions and actuator power required are performed 
for the cases when only individual components of the vibratory hub shears are suppressed. For this 
study, the fixed-gain controller is used to suppress the vibratory hub shears for the same three-bladed 
articulated rotor considered previously. The flight speed is held at 135 knots at a thrust coefficient- 
solidity ratio of 0.08. Figure 3.9 shows the comparison between the controlled vibratory hub shears, 
when each component is suppressed individually. For comparison purposes, the uncontrolled values 
are also presented. The comparisons of HHC swashplate amplitudes and actuator powers are shown in 
figure 3.10. The inputs required to suppress all three vibratory hub shear components are also shown. 
The suppression of the 3 P vertical hub shear is achieved with only the collective input, corresponding 
to the third harmonics of the blade pitch schedule. When this vibratory hub shear component is 
completely suppressed, there are small increases in the 2 P and the 4 P rotating inplane hub shears from 
the uncontrolled values. Although the collective swashplate amplitudes required for this case are nearly 
equal to those required for the complete suppression of all three vibratory hub shears, the actuator power 
required is much smaller. 

For the case when the 2 P rotating inplane hub shear is suppressed, only the cyclic swashplate 
motion corresponding to the second harmonic blade pitch is employed. Figure 3.9 shows that there 
is a small increase in the 3 P vertical hub shear and a small decrease in the 4 P rotating inplane hub 
shear due to this HHC component. It is worth noting from figure 3.10 that suppression of the 2 P 
rotating inplane hub shear requires approximately the same amount of cyclic actuator input, as well as 
the actuator power, as that required for the suppression of all three vibratory hub shear components. 
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Finally, suppression of the 4 P rotating inplane hub shear using only the cyclic swashplate motion, 
corresponding to the fourth harmonics of blade pitch, is accompanied by about 30 percent reduction in 
the 3 P vertical hub shear and a small decrease in the 2 P rotating inplane component (fig. 3.9). The 
swashplate amplitudes and the actuator power required to suppress this vibratory hub shear component 
are slightly smaller than those required to suppress only the 3 P vertical hub shear. 

Even though the baseline value for the 3 P vertical hub shear is much larger than the 2 P rotating 
inplane component (fig. 3.9), the actuator power required to suppress the former component is about 
one-quarter of that required for the suppression of the latter component. This is understandable since 
HHC uses the higher harmonic blade pitch to generate new airloads to suppress the vibratory hub loads. 
Since a rotor blade generates more lift (in the vertical direction) than drag (in the inplane direction) for 
a given change in blade angle of attack, the vertical hub shear is more sensitive to the higher harmonic 
inputs than the inplane components. Therefore, larger HHC inputs are required to suppress the 2 P 
rotating inplane hub shear. 

The validation study for the three-bladed articulated model rotor of the CH-47D is also performed 
for cases when the rotor thrust is varied. For this, the forward speed is held fixed at 60 knots corre- 
sponding to an advance ratio of 0.138. The rotor thrust coefficient-solidity ratio (C T /a) is varied from 
0.04 to 0.14. This flight condition generally represents the transition flight regime. Figure 3.11 shows 
the correlation of the 3 P vertical hub shear with and without HHC effects. With no HHC and at a low 
thrust level (C T / a = 0.04), the analytical results agree quite well with experimental data. For values of 
Or / ,<J greater than 0.08, the calculated value of the 3 P vertical hub shear overpredicts the experimental 
data. This vibratory hub shear component increases with thrust level, and the trend is captured well 
in the analytical results. The discrepancies at high thrust may be attributed to the differences between 
the calculated trim controls and the experimental values. These differences are small for the collective 
controls, whereas the differences in cyclic controls are as large as 3 degrees. 

Figure 3. 12 shows the correlation of the 2 P rotating inplane hub shear for the same flight conditions. 
For most of the range of thrust level shown, the analytical results underpredict the experimental data. 
However, the increasing trend of the vibratory hub shear component with increasing rotor thrust is well 
predicted with the analysis. Correlation of the 4 P rotating inplane hub shear is shown in figure 3.13 and 
indicates that the calculations underpredict the experimental data for the complete range of thrust. The 
analytical results indicate that this vibratory hub shear component, being much smaller as compared to 
other components, increases slightly with thrust level. It was also found that perturbations in the cyclic 
trim values have a substantial influence on the predicted vibratory inplane hub shears. The discrepancies 
between the predicted and the measured 4 P rotating inplane hub shears may be attributed to the rotor 
inflow behavior at this flight condition, especially for wind tunnel testing. For a rotor tested in a wind 
tunnel at the transition flight regime, a significant source of unsteady airloads originates from the rotor 
wake which remains close to the rotor disk. Furthermore, this problem is exacerbated because the 
rotor wake geometry is also influenced by the presence of the wind tunnel walls and floor, and such 
recirculation effects have not been considered in the simulation study. Another factor which might affect 
the results of the correlation is the effect of the hub motion which was not included in the analysis. 
Compared to the hub loads predicted by other authors, the present analytical results are in much better 
agreement with the experimental values. 
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For the same flight conditions, HHC was introduced into the analysis to suppress the predicted 
vibratory hub shears. The fixed-gain HHC controller with no control rate limiting factor (r = 0) was 
used. Again, both analytical and experimental results are presented in figures 3.11-3.13. These figures 
show that the fixed-gain controller suppresses almost completely all three vibratory hub shear compo- 
nents. This effectiveness is confirmed by both analytical and experimental results. Since measurement 
and process noise were not considered in the analysis, the predicted vibratory hub shears are suppressed 
to a lower level than those obtained experimentally. 

The corresponding swashplate amplitude and phase variations with thrust level are shown in 
figures 3.14 and 3.15, respectively. These results indicate that all components of the HHC swashplate 
amplitudes increase with the thrust level. The variation is almost linear for the collective swashplate 
amplitude. It is interesting to note that the phases of both cyclic actuators vary by about 60 degrees as 
Cx/cr is varied from 0.04 to 0.14 while the collective actuator phase remains relatively constant. The 
actuator power required by the HHC system, expressed in terms of the actuator power index, is shown 
in figure 3.16. Within the range of rotor thrust shown, the actuator power increases quite significantly 
with the rotor thrust level. The actuator power required at the high thrust level ( Op/cr = 0.14) is about 
four times that for the low thrust level (C T /cr = 0.04). 

3.1.2 Secondary Effects of HHC 


Next, the effects of HHC on the alternating blade bending stresses, control system loads, and the 
rotor shaft torque are examined (figs. 3.17—3.21). Alternating blade bending stresses are quite important 
since these affect the structural life of rotor blades. In figures 3.17—3.19, the alternating blade bending 
moments near the blade midspan for the flap, lag, and torsion modes are shown, respectively, with and 
without HHC effects. Due to the applied HHC, there are reductions in both the flap and the lag bending 
moments at 0.48 and 0.31 radius, respectively (figs. 3.17 and 3.18). Interestingly, the magnitudes 
decrease with thrust level and can be as much as 30 percent less than the baseline values. 

In contrast, the alternating blade torsional moment at the blade midspan increases during HHC 
application over the range of thrust levels shown in figure 3.19. The increase in this blade stress 
component can be as much as 50 percent of the baseline values. Figure 3.20 shows the alternating pitch- 
link load as a function of thrust level. With HHC, it increases up to 48 percent from the uncontrolled 
values. An increase in the alternating blade torsional moment at midspan and the control loads were 
also observed in wind tunnel testing. These increases are expected since the HHC inputs are in fact 
implemented by exciting the blade torsional modes. 

The effects of HHC on rotor performance, presented in terms of the rotor torque in figure 3.21, 
are negligible. Changes in the rotor shaft torque with HHC on range from a 0.35 percent increase at 
the low thrust level ( C^/cr = 0.04) to a 0.5 percent reduction at the high thrust level (Cj’/cr = 0.14). 
Negligible changes in the rotor shaft torque with HHC were also noted in the wind tunnel investigation 
(ref. 46). 
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3.1.3 HHC Induced Blade Stall and Controller Performance 


In order to investigate the influence of HHC on blade stall, a flight condition at 135 knots ( (i = 0.3) 
and Cj'/a — 0.09 is used. For this, the global (fixed-gain) controller (r = 0) is used to suppress the 
three vibratory hub shear components: 3 P vertical and the 2 P and the 4 P rotating inplane hub shears. 
The vibratory responses presented in figure 3.22 show that all three components of the vibratory hub 
shear are completely suppressed within three HHC cycles. As an indicator of blade stall, the blade 
section normal force coefficients computed at 93 percent blade radius for the cases with HHC off and 
HHC on are shown in figure 3.23 and 3.24, respectively. In these figures, C' N is the stall monitor 
value of the normal force coefficient, Cpj is the total normal force coefficient (including dynamic stall 
effects), and the stall boundary is the critical normal force value, Cjy l . As discussed in the formulation, 
dynamic stall is initiated when the value of C' N exceeds the critical normal force value, C N . A stall 
margin is thus defined as the difference between the values of C' N and C]\f . 

Figure 3.23 shows that, with no HHC, the stall margin value for this particular blade section near 
the blade tip is roughly 0.4 to 0.5 for the complete range of the rotor azimuth. However, figure 3.24 
shows that, with HHC on, the stall margin for this particular blade section is reduced, especially in 
the retreating side of the disk in the fourth quadrant. This result, although limited to a particular blade 
section, indicates that the HHC inputs used to alleviate helicopter vibration might promote blade stall. 

Further investigation of the HHC effects on blade stall is made by raising the Cj'/a to 0.135 for 
the same forward speed. This hypothetical flight condition exceeds the level-flight structural envelope 
and maneuver capability of the CH-47 helicopter (ref. 99). Again, the global controller (r = 0) is first 
used to suppress the vibration, and the resulting vibratory hub shear responses are shown in figure 3.25. 
Clearly, complete vibration suppression for this flight condition can not be achieved using the global 
controller. Over the range of HHC cycles, the vibratory hub shears have essentially the same magnitudes 
as the uncontrolled values. 

In figures 3.26 and 3.27, the resulting section normal force coefficients near the blade tip (93 percent 
radius) for the case discussed above are shown without and with HHC over one rotor revolution. The 
results presented in figure 3.27 are obtained at the sixth HHC cycle. Without HHC, the blade section 
is already stalled in the fourth quadrant, whereas, with HHC on, there are two large lift overshoots that 
cross over the stall boundary in the retreating blade tip region (fig. 3.27). These results indicate that for 
a rotor operating with considerable stall on the rotor disk, the HHC inputs used to suppress the vibratory 
hub loads escalate the stall phenomena and thus deteriorate the HHC performance. 

For the same flight condition, the local model with a control rate limiting factor, r = 0.3, is then 
used to suppress the vibration. Figure 3.28 shows the time history of the vibratory hub shear responses. 
This figure shows that the local controller, although having a sluggish transient behavior (within the 
first two cycles), is quite effective in suppressing all three vibratory hub shear components. At the 
eighth HHC cycle, there is a reduction of more than 95 percent for all three hub shear components. 
This controller effectiveness is attributed to two effects. First, the transfer matrix at this flight condition 
is nonlinear, and the local model is more suitable for such conditions. Second, the control rates are 
constrained to suppress large control excursions, and this is beneficial for avoiding excessive blade stall. 
The resulting blade section normal force coefficients at 93 percent blade radius for the eighth HHC 


66 



cycle are shown in figure 3.29. These results should be compared with those presented in figure 3.26. 
Even though this controller promotes some blade stall on the retreating side of the disk, the resulting 
lift variation appears much smoother than that obtained with the global controller. 

For the same forward speed, the rotor thrust is now reduced to the C T /cr = 0.12 level. This flight 
condition represents the level-flight boundary condition of the CH-47 aircraft (ref. 99). For this flight 
condition, the vibration is suppressed within four cycles using the local controller with a control rate 
limiting factor of 0.3. The effects of HHC in promoting blade stall, viewed from above the rotor disk, 
are shown in figures 3.30 and 3.31. In these figures, the stall area on the disk is hatched and corresponds 
to the locations where the value of C' N is greater than the C Nl value. Figure 3.30 shows that with no 
the rotor is already stalled in the third and fourth quadrants of the disk, as might be expected. 
However, the retreating blade tip region is stall free. Figure 3.31 indicates that, with HHC on, the stall 
area spreads up to the blade tip region in the fourth quadrant of the disk. Again, this result confirms 
that HHC inputs used to suppress vibration can indeed promote blade stall. 

As a final note in this investigation, the results presented are limited to the wind tunnel simulation 
of the three-bladed articulated rotor model of the CH-47D helicopter. For other rotors with different 
blade numbers or with a different hub design, such as hingeless or bearingless, the results may be 
entirely different. The HHC inputs used to suppress the rotor vibratory hub loads might not induce 
blade stall on the rotor disk. Furthermore, the trim conditions used in this study simulate the flight 
conditions of the CH-47D which is a tandem helicopter. For a single main rotor helicopter, the trim 
requirements are entirely different, and a different HHC-stall-induced behavior may be expected. 


3.2 Sensitivity Analysis 


3.2.1 Hingeless Versus Articulated Rotors 


The above investigation was performed for a typical articulated rotor system. Next, the vibration 
suppression characteristics of two hingeless rotor configurations are analyzed and compared with results 
of the articulated rotor system. These three-bladed hingeless rotors represent, respectively, a soft-inplane 
(lag frequency = 0.7 IP) and a stiff-inplane (lag frequency = 1.4P) configuration. Both rotors have a 
fundamental flap frequency of 1 . 1 P. For these results, wind tunnel flight conditions are simulated with 
the forward speed ranging from zero to 188 knots and at a thrust coefficient-solidity ratio of 0.08. 

Figure 3.32 shows the variation of the 3P vertical hub shear with forward speed for three rotor 
systems. Below 100 knots, the calculated 3 P vertical hub shears for these rotors are quite comparable. 
At the higher forward speeds, the 3 P vertical hub shear for the stiff-inplane rotor increases rapidly with 
speed, while this component remains the lowest for the articulated rotor system. 

Variation of the 2 P rotating inplane hub shear with forward speeds for these three rotors is shown in 
figure 3.33. For forward speeds below 135 knots, the 2 P rotating inplane hub shears for the articulated 
and the soft-inplane hingeless rotors are quite similar, whereas at higher speeds, the vibratory hub shear 
component for the soft-inplane hingeless rotor is about 30 percent lower than that for the articulated 
rotor system. For the stiff-inplane rotor, the 2 P vibratory hub shear increases considerably for forward 
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speeds above 100 knots. This is understandable since the blade fundamental lag mode frequency for the 
stiff-inplane rotor of 1.4P is closer to the resonant frequency of 2 P than the other two rotor systems. 
Figure 3.34 shows the 4 P rotating inplane hub shears. The magnitudes of this high frequency hub shear 

are much smaller than those of the low frequency 2 P hub shear and are quite similar for ail three rotor 
systems. 


Next, HHC inputs are applied to suppress the vibratory hub shear components for these three rotor 
configurations. Because the HHC system is quite effective in suppressing completely the 3 P vertical 
and the 2 P and 4P rotating inplane hub shears, the controlled vibratory hub shears are not shown. Even 
though the vibratory hub moments for the hingeless rotor are hot small, especially at the high forward 
speeds, these are not suppressed in this study, so subsequent HHC results can be compared directly with 
those obtained earlier for the articulated rotor. 


Figure 3.35 shows the actuator power index required by the HHC system at different forward 
speeds. For the stiff-inplane hingeless rotor above 100 knots, the HHC inputs required to suppress 
the vibratory hub shears become very large, and the problem associated with the convergence of the 
numerical solution hinders any more results. For moderate to high forward speeds, the power index for 
the soft-inplane hingeless rotor is larger than that required for the articulated rotor. It is interesting to 
note that in the transition flight regime, the actuator power index is the largest for the soft-inplane rotor 
even though the vibratory hub loads are the largest for the stiff-inplane configuration. However, for 
forward speeds above 140 knots, the power index for the soft-inplane hingeless rotor becomes lower 
than that for the articulated rotor system. In fact, at 160 knots, the power index for the soft-inplane 
hingeless rotor is less than half that of the articulated rotor. For forward speeds above 50 knots, the 
actuator power index required by the stiff-inplane hingeless rotor is the largest among the three rotor 
systems considered. 

3.2.2 Parametric Study— Four-BIaded Hingeless Rotor 


Next, the sensitivity of the uncontrolled vibratory hub loads and the actuator power for a hingeless 
rotor with respect to various blade design parameters is examined. For this simulation, a four-bladed 
soft-inplane hingeless rotor is chosen as a baseline configuration. The structural characteristics of this 
rotor blade are similar to those of the three-bladed soft-inplane hingeless rotor considered earlier, i.e., 
the blade fundamental flap and lag mode frequencies are 1.1P and 0.7P, respectively. Propulsive trim 
with steady level flight is simulated for forward speed of 40 knots corresponding to a rotor advance 
ratio of 0.09. This airspeed corresponds to the transition flight regime where the rotor would experience 
a high level of vibration. The rotor thrust coefficient-solidity ratio is held at 0.08. The fixed-gain HHC 
controller is now configured to suppress the 4P vertical, the 3 P and 5 P rotating inplane hub shears, 
and the 3 P and 5P rotating hub moments. 

Blade torsion stiffness- The blade torsion frequency is varied by adjusting the blade torsion 
stiffness of the inboard section. In this simulation, blade rotating fundamental torsion frequency is 
varied from 3 P to 7 P, and the baseline value is 4.5 P. Figure 3.36 shows the uncontrolled vibratory 
hub shears with blade torsion frequency. For most of the range of frequencies shown, the rotor vibratory 
hub shear components are quite insensitive to the blade torsion stiffness. An exception is for frequencies 
close to 3 P where there is a considerable reduction in the 4 P vertical hub shear and a moderate increase 
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in the 3 P inplane hub shear. This is attributed to the resonance of blade torsional response with the 
forcing at 3 P. The vibratory hub moments are presented in figure 3.37, which shows that the 3 P hub 
moment increases moderately with blade torsion stiffness. The 5 P hub moment is small compared to 
the 3 P component and is quite insensitive to blade torsion stiffness. The variation of the actuator power 
index with blade torsion frequency is shown in figure 3.38. It is interesting to note that the power index 
is lowest at the blade torsion frequency of about 3.85P, for which the actuator power index is less than 
half that at the torsion frequency of 7 P. 


Blade bending stiffness- Blade bending stiffnesses in the flapwise and chordwise directions are 
varied to study their effects on the vibratory hub loads and the actuator power. Blade bending stiffness 
in the flapwise direction affects the blade elastic flap modes, and its variation is presented in terms of 
the blade second rotating flap mode frequency. Figure 3.39 shows the uncontrolled vibratory hub shear 
variation with blade flap mode frequency; the baseline value is 2.85P. Variation of the flapwise stiffness 
affects all three vibratory hub shear components. The 4 P vertical hub shear is a minimum when the 
blade frequency is about 3 P , while the 3 P inplane hub shear is a minimum when the blade frequency 
is 2.65 P. The 5 P inplane hub shear component is small compared to the other two components, and 
its magnitude decreases steadily with an increase in flap stiffness. Variations of the 3 P and 5 P hub 
moment with blade flap bending stiffness are shown in figure 3.40. Within the range of flap frequencies 
shown, it is interesting to note that the 5 P hub moment increases sharply as the flap frequency is reduced 
below 2.6P, while the 3 P hub moment reaches a minimum at this frequency value. In figure 3.41, the 
actuator power remains constant for the range of flap frequencies above 2.9 P and increases sharply for 
the range of flap frequencies below this value. 

Variation of the blade chordwise stiffness is presented in terms of the blade second rotating lag 
mode frequency; the baseline value is 5.8P. The uncontrolled vibratory hub shears are presented in 
figure 3.42. For this rotor configuration, the vibratory hub shears are quite sensitive to the placement of 
the blade second lag mode frequency. The magnitude amplification of hub shears is quite pronounced 
as the blade chordwise frequency is placed at the corresponding harmonic. It is interesting to note that 
the resonant conditions exist not only for the inplane shears but also for the vertical component. In fact, 
at resonance, the 4P vertical hub shear increases by four and a half times the baseline value. Magnitude 
amplification is also shown in figure 3.43 where the vibratory hub moments are presented. This figure 
shows that the 5P hub moment increases noticeably as the lag frequency approaches 4 P, while the 
effect of lag stiffness on the 3P hub moment is small. The actuator power index variation shown in 
figure 3.44 indicates that minimum power required is achieved by placing blade lag frequencies either 
at 3.6P or at 6.2 P away from the resonant frequencies at the rotor harmonics. 

Blade linear pretwist- Blade pretwist redistributes blade steady airloads in the spanwise direction 
and introduces elastic flap-lag coupling. In this study, blade linear pretwist is used with values ranging 
from zero to -16 degrees. The baseline blade configuration has a pretwist of -8 degrees. Figures 3.45 
and 3.46 show, respectively, the vibratory hub shears and moments with blade pretwist. These figures 
show that the 4P vertical and the 3 P inplane hub shears and the two hub moment components all 
decrease with increasing blade pretwist. The 5P inplane hub shear is small and increases slightly with 
increasing blade pretwist. The actuator power is shown in figure 3.47, where a reduction of 20 percent 
in the actuator power is achieved with a blade pretwist of -16 degrees as compared to an untwisted 
blade. 
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Blade chordwise center-of-gravity offset- Blade chordwise center-of-gravity (c.g.) offset from 
the elastic axis introduces inertial coupling between the flap and torsion modes. For this study, the elastic 
axis is assumed to be at the quarter chord, and a positive c.g. offset represents the blade center of mass 
placed ahead of the quarter chord, toward the blade leading edge. The baseline rotor configuration has 
no offset. Figures 3.48 and 3.49 show that the 4 P vertical and the 3 P inplane hub shears and the 3P 
hub moment are all reduced as the blade c.g. is shifted forward of the elastic axis. The 5 P inplane hub 
shear and the 5 P inplane hub moment are small compared to the other vibratory hub load components 
and increase slightly as the blade chordwise c.g. is moved ahead of the elastic axis. From figure 3.50, 
the actuator power increases sharply as the blade c.g. is placed more than 2 percent of chord ahead of 
the elastic axis and increases slightly as the blade c.g. is placed behind the elastic axis. In fact, the 
actuator power is lowest when the blade c.g. is placed at 2 percent of chord length, ahead of the elastic 
axis. 


Offset of elastic axis from quarter-chord axis- The offset of blade elastic axis from the blade 
quarter chord introduces aerodynamic pitch-flap coupling, and its effects on the vibratory hub loads and 
HHC system were examined. The offset is positive when the blade elastic axis is placed ahead of the 
quarter chord, and, for the baseline rotor, the elastic axis is at quarter chord (zero offset). Figures 3.51 
and 3.52 show the variations of the vibratory hub loads with the variation of this offset value. Variations 
of the vibratory hub loads are quite similar to those observed when the blade c.g. is varied with respect 
to the elastic axis location. Figure 3.53 shows the actuator power variation with location of elastic axis, 
and a similar behavior is observed to that shown in figure 3.50. Placing the elastic axis in front of the 
quarter chord helps reduce the major vibratory hub load components. A minimum in actuator power 
is achieved by placing the elastic axis 3 percent of chord length ahead of the quarter chord. At this 
minimum value, a 17 percent reduction in actuator power is achieved from the baseline value. 

Blade thrust- Rotor thrust, expressed in terms of C T /a, is varied from 0.06 to 0.12, and the results 
are shown in figures 3.54—3.56. Figures 3.54 and 3.55 show that all the vibratory hub load components 
increase noticeably with an increase in rotor thrust. The actuator power variation is presented in 
figure 3.56, which shows that, at high thrust {Cj'fa = 0.12), the power required is twice that for the 
baseline case (Cj 1 / cr = 0.08). Within the range of rotor thrusts shown, the actuator power varies almost 
linearly with the rotor thrust. 

3.2.3 Parametric Study— Three-Bladed Articulated Rotor 


Next, the sensitivity of the uncontrolled vibratory hub loads and the actuator power for the baseline 
three-bladed rotor with respect to various blade design parameters is examined. For this simulation, the 
level-flight boundary condition of the aircraft is used, i.e., V = 135 knots (/x = 0.3), and C T /a - 0.12. 
The closed-loop HHC controller with control rate limiting factor is used to suppress the 3 P vertical and 
the 2 P and 4 P rotating inplane hub shears. The blade and rotor properties considered previously in the 
correlation study (table 3.1) are used as a baseline rotor configuration. 

Blade torsion stiffness— In this simulation, the blade fundamental torsion frequency is varied from 
3.8P to 6.8P; the baseline value is 4.86P. The uncontrolled vibratory hub shear variation with blade 
torsion frequency is presented in figure 3.57. For most of the frequency range shown, the vibratory 
inplane hub shears are quite insensitive to the blade torsion frequency. For the torsion frequencies below 


70 



the baseline value, the 2 P rotating inplane hub shear increases slightly, while the 4 P inplane component 
decreases slightly from the baseline values. The 3 P vertical hub shear increases moderately with blade 
torsion stiffness. In fact, there is a 35 percent reduction in the 3 P vertical component as the blade 
torsion frequency is reduced from 6.8P to 3.8P. The actuator power index variation with blade torsion 
frequency is shown in figure 3.58. From an actuator design point of view, it is interesting to note that a 
torsionally stiff blade requires less actuator power than a torsionally soft blade. A torsionally stiff blade 
with frequency greater than 4.5P consumes roughly 40 percent less actuator power than a blade with 
torsion frequency placed at 3.8P. For torsion frequencies greater than 4.5P, the actuator power index 
is insensitive to the blade torsion vibration characteristics. 

Blade bending stiffnesses- Blade bending stiffnesses in the flapwise and chordwise directions are 
varied to study their effects on the vibratory hub loads and the actuator power. Figure 3.59 shows 
the uncontrolled vibratory hub shear variation with the blade flap mode frequency, where the baseline 
frequency value is 2.66P. Variation of flapwise stiffness significantly affects the 3 P vertical hub shear, 
while the inplane components are less sensitive to this parameter. The value of the 3P vertical hub shear 
is highest for low flap stiffness. As the flap frequency is increased, this vibratory hub shear component 
reaches a minimum at 2.78 P, increases to a local maximum at 3 P, and then decreases again for flap 
frequencies greater than 3P . The peak value which occurs at 3 P may be attributed to the resonant 
condition. The 2 P rotating inplane hub shear increases slightly with increasing flap stiffness, while the 
4P component remains relatively unchanged. The actuator power index variation with flap frequency 
is shown in figure 3.60. This figure indicates that minimum actuator power is achieved at the first flap 
bending frequency of 2.65P . Reducing the flap stiffness from this optimum value does not significantly 
change the actuator power required, whereas increasing it would increase the actuator power by a small 
amount. 

Variation of the blade chordwise stiffness is presented in terms of the blade first lag bending 
frequency, the baseline frequency is 6.4P. The uncontrolled vibratory hub shears are presented in 
figure 3.61. The sharp increase in the 4P rotating inplane hub shear as the blade lag bending frequency 
approaches 4 P is attributed to the resonant condition. At this frequency, the 4P inplane hub shear 
increases by five times the baseline value. For the range of frequencies greater than 4.5P, this vibratory 
hub shear component is quite insensitive to the blade lag stiffness. Within the range of blade lag 
frequencies shown, the 2 P rotating inplane hub shear is insensitive to blade lag stiffness, while the 3P 
vertical component varies significantly. At the lag frequency of 4P, the 3 P vertical hub shear is highest. 
As the lag stiffness is increased, the 3 P vertical hub shear component reaches a local minimum at 4.6P, 
increases to a local maximum at 5.7P, decreases again to another local minimum at 6.4P, and again 
increases at the higher frequency values. The actuator power index variation is shown in figure 3.62. 
The HHC power required to suppress vibration is quite insensitive to the blade lag stiffness, except 
when the lag frequency is placed at 4P, at which point the actuator power is increased by 30 percent 
from the baseline value. 

Blade chordwise c.g. offset- The baseline rotor configuration has a nonzero c.g. offset distribution 
as shown in table 3.1. Figure 3.63 shows the uncontrolled vibratory hub shear as a function of blade 
chordwise c.g. offset from the elastic axis, expressed in percentage of chord length. Except for the 
baseline configuration, the chordwise c.g. offset considered in this parametric investigation is made by a 
uniform variation across the blade span. Figure 3.64 shows that all three vibratory hub shear components 
are reduced as the blade c.g. is moved ahead of the elastic axis. The most significant reduction is the 
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2P rotating inplane hub shear which is reduced by more than half as the blade c.g. is shifted uniformly 
by more than 6 percent of chord length ahead of the elastic axis. Figure 3.64 shows that the actuator 
power is reduced as the blade chordwise c.g. is placed ahead of the elastic axis. In fact, by moving the 
blade c.g. from the elastic axis toward the blade leading edge by 10 percent of chord length reduces the 
actuator power index by 40 percent. 

Offset of elastic axis from quarter chord- Figure 3.65 shows the uncontrolled vibratory hub 
shear variation with the offset of elastic axis from quarter chord. The 3 P vertical and the 4 P inplane 
hub shears are quite insensitive to this parameter, while the 2 P inplane hub shear varies significantly 
with the elastic axis location. As the elastic axis is shifted forward of the blade quarter chord, the 2 P 
inplane component reduces to a minimum at 5 percent of chord length offset and increases at higher 
offset value. Figure 3.66 shows the actuator power index variation with the offset of the elastic axis 
from quarter chord. From this figure, the actuator power is lowest when the elastic axis is placed at 
5 percent of chord length ahead of the quarter chord, at which point the actuator power is reduced by 
half from the baseline value. 
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4 CONCLUSIONS 


An analytical formulation has been developed to calculate the vibratory hub loads for helicopter 
rotors in forward flight and the corresponding higher harmonic control inputs to suppress these loads. 
The analysis was based on a finite element method in both space and time domains. A nonlinear time 
domain unsteady aerodynamic model, based on the indicial response formulation, was used to calculate 
the blade airloads. The rotor induced inflow was calculated using a free wake analysis. The vehicle 
trim control settings and blade steady responses were solved as one coupled solution using a modified 
Newton method. The HHC responses form a part of the steady response calculations where the basic 
collective and cyclic trim inputs were superimposed with the higher harmonic components. A linear 
frequency-domain quasi-steady model, in terms of a transfer function, was used to relate the harmonics 
of the HHC inputs to the harmonics of the vibratory hub loads. Optimal HHC was calculated from the 
minimization of the vibratory hub loads, expressed in terms of a quadratic performance index. 

The calculated vibratory hub loads, with and without HHC effects, have been correlated with wind 
tunnel data obtained from a three-bladed articulated scaled model rotor. The secondary effects of HHC 
were also investigated including the influences of HHC on blade stall. Systematic parametric studies 
have been performed on different rotor systems to investigate the effects of blade design variables on 
the vibratory hub loads, HHC inputs, and actuator power requirements. The following conclusions were 
drawn from this study: 

1 . Good correlation of predicted vibratory hub loads was obtained with experimental data at 
different flight conditions. A refined aerodynamic model including compressibility effects and dynamic 
stall helped improve the correlation at high forward speeds. The inclusion of the free wake analysis 
improved the correlation in the transition flight regime. For the three-blade model rotor operated at low 
speed, the correlation of the 3 P vertical and the 2 P rotating inplane hub shears was adequate, while 
the calculated 4 P rotating inplane hub shear was significantly lower than the experimental values. 

2. The fixed-gain HHC controller was quite effective in suppressing the rotor vibratory hub shears 
for most of the rotor operating conditions. The HHC actuator amplitudes and power requirements 
increased substantially at high forward speed range (above 100 knots, or fi = 0.2) and increased 
moderately with the rotor thrust level. For complete vibratory hub shear suppression, a large percentage 
of the actuator power is required to suppress the 2 P rotating inplane hub shear component on a three- 
bladed rotor. At high forward speeds, a soft-inplane hingeless rotor required less actuator power than an 
articulated rotor, whereas the actuator power requirement for a stiff-inplane hingeless rotor was generally 
higher than that for an articulated or a soft-inplane hingeless rotor. 

3. For a rotor operating in the transition flight regime, HHC caused moderate reduction in the flap 
and lag bending stresses near the blade midspan, whereas it increased the control system fatigue loads 
and the blade torsional stresses significantly. HHC effects on overall rotor performance were found to 
be negligible at the low speed regime. 

4. For a rotor operated at or outside the flight boundary, the adaptive-gain controller suppressed 
almost completely the rotor vibratory hub shears, even though HHC might actually increase the stalled 
flow region on the rotor disk. For such flight conditions, the fixed-gain controller was ineffective in 
suppressing the vibratory hub loads. 
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5. For a four-blade soft-inplane hingeless rotor operated in the transition flight regime 
(/X = 0.09, C T /cr = 0.08), 

a. The actuator power was quite sensitive to the blade torsion stiffness, while the baseline 
vibratory hub loads were less sensitive. Placing blade torsion frequency at 3.85P resulted 
in a minimum actuator power. 

b. Placing blade second elastic flap frequency above 2.9P helped reduce the actuator power. 

c. The uncontrolled vibratory hub loads were sensitive to the placement of the blade second 
elastic lag frequency. Placing blade second elastic lag frequency at 3.6P or at 6.2P helped 
reduce the actuator power. 

d. An increase in blade linear pretwist reduced the actuator power. 

e. Placing blade center of mass 2 percent of chord length ahead of the elastic axis and placing 
the elastic axis 3 percent of chord length forward of the quarter-chord axis helped reduce 
the actuator power. 

f. The uncontrolled vibratory hub loads and the actuator power increased with rotor thrust. 

6. From a sensitivity study of blade design parameters for a three-bladed articulated rotor operating 
at the level-flight boundary (/x = 0.3, C T /a = 0.12), 

a. The 3 P vertical hub shear was sensitive to the blade torsion stiffness, while the inplane 
components were less sensitive. The actuator power was quite insensitive to the blade 
torsion frequency, except near the resonant frequency of 3 P. 

b. The 3 P vertical hub shear was sensitive to the placement of the blade first flap bending 
frequency, while the inplane vibratory hub shear components and the actuator power were 
less sensitive. Both the 3 P vertical and 4 P rotating inplane hub shear were sensitive to 
the blade lag stiffness. Placing the blade first elastic lag frequency at AP increased the 
AP rotating inplane hub shear by five times and the actuator power by one-third from the 
baseline values. 

c. Placing the blade center of mass ahead of the elastic axis helped reduce both the uncontrolled 
vibratory hub shears and the actuator power. Placing the elastic axis at 5 percent of chord 
length ahead of the quarter chord minimized the 2 P rotating inplane hub shear as well as 
the actuator power. 


Items (5) and (6) above indicate that blade parameters had moderate to strong influence on the 
amplitude and actuator power requirements for an HHC system. Therefore, the design of an optimum 
blade which is best suited for HHC application (i.e., with minimum actuator power requirement) is, in 
principle, feasible. 
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Table 3.1. Rotor and blade properties 


Element 

Ely 

EI Z 

GJ 

m 

m 0 

4 . 10 3 

R 2 

EL 

c 

/ 

m 0 ti 2 R 4 

mo0 2 7? 4 

m 0 Q 2 R 4 

72 

1 

0.002596 

0.07536 

0.001331 

1.0000 

0.3312 

0.03226 

0.1500 

2 

0.003326 

0.07782 

0.002538 

0.8199 

0.2779 

0.02268 

0.1500 

3 

0.003744 

0.08020 

0.002546 

0.8807 

0.2728 

0.03194 

0.2000 

4 

0.003966 

0.08020 

0.002578 

1.2065 

0.2528 

0.04616 

0.2200 

5 

0.006306 

0.05961 

0.003470 

0.9938 

0.1269 

0.01694 

0.1352 

6 

0.016667 

0.01785 

0.014480 

4.5808 

0.2433 

0.00385 

0.1162 


Blade number, N 3 

Blade radius, 72 5.23 ft 

Rotor rotational speed, Cl 22.226 Hz 

Reference blade mass per unit length, mo 0.0322 lb/in. 

Chord length/radius, c/R 0.07446 

Solidity ratio, a 0.0711 

Pitch-link stiffness, Kp/m^CfiR 3 0.015165 

Flap hinge offset 0.028617 

Lag hinge offset 0.14487? 

Root cut-out 0.2172 

Blade linear pretwist, 9tw —14 deg 

Lock number, 7 = SpcRqJmo 6.8468 

Lift curve slope, c ia at M = 0.3, per radian 6.1306 

Viscous drag coefficient, c ^ 0 0.0085 

Leading edge suction recovery factor, 77 0.97 

Pitching moment coefficient at zero angle of attack, c mo -0.005 
Hover tip Mach number 0.6415 


Table 3.2. Blade natural frequencies 


Mode 

Blade natural 
frequency, per rev 

Rigid lag 

0.504 

Rigid flap 

1.023 

First elastic flap 

2.660 

Rigid torsion 

4.865 

Second elastic flap 

5.210 

First elastic lag 

6.391 

Third elastic flap 

8.918 

First elastic torsion 

11.439 
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Blade spanwise lift distribution, 
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— Analytical, Nonlinear Aero., Free Wake 





Figure 3.2. Correlation of 2 P rotating inplane hub shear — variation with forward speed (C^/cr = 0.08). 



Figure 3.3. Correlation of 4 P rotating inplane hub shear — variation with forward speed (Cj'/cr — 0.08). 
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NORMALIZED 4P ROTATING 
INPLANE HUB SHEAR 



Figure 3.6. HHC suppression of 4 P rotating inplane hub shear ( C T /a = 0.08). 



Figure 3.7. Variation of HHC swashplate amplitudes with forward speed (Cj>/g = 0.08). 
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Figure 3.8. Variation of actuator power index with forward speed (C T /a ~ 0.08). 
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Figure 3.9. Effects of closed-loop HHC controller on the suppression of individual vibratory hub shear 
components ( V = 135 knots, C T /a =0.08). 
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Figure 3.10. Actuator amplitudes and power required for the suppression of individual components of 
vibratory hub shear (V = 135 knots, C T /a = 0.08). 



Figure 3.11. Correlation of 3 P vertical hub shear with and without HHC effects — variation with rotor 
thrust (V = 60 knots). 
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Figure 3.12. Correlation of 2 P rotating inplane hub shear with and without HHC effects — variation with 
rotor thrust (V = 60 knots). 
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Figure 3.13. Correlation of AP rotating inplane hub shear with and without HHC effects— variation with 
rotor thrust (V = 60 knots). 
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Figure 3.14. Variation of actuator amplitudes with thrust level ( V = 60 knots). 
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Figure 3.15. Variation of actuator phase angle with thrust level (V = 60 knots). 
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Figure 3.16. Variation of actuator power index with thrust level ( V = 60 knots). 
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Figure 3.17. HHC effects on the alternating flap bending moment (V = 60 knots). 
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Figure 3.18. HHC effects on the alternating lag bending moment (V = 60 knots). 


om 

</> y 

o 5 30 

I— w 

U1QC 

QO 

~i • 20 

CD . 


HHC Off 

HHC On 


V » 60 Knots 
(M = 0.14) 


< LU IQ 
22 
DC O 


O' 


THRUST COEFFICIENT-SOLIDITY RATIO, C t /ct 

Figure 3.19. HHC effects on the alternating torsional moment (V = 60 knots). 
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Figure 3.20. HHC effects on the alternating pitch-link load (V = 60 knots). 
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Figure 3.21. HHC effects on the rotor shaft torque (V — 60 knots). 
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Figure 3.22. Vibratory hub shear responses with fixed-gain HHC controller, r = 0.0 (V = 135 knots 
C T /a = 0.09). 
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Figure 3.23. Variation of normal force coefficient at 0.93 radius with HHC off ( V = 135 knots 
C T /a = 0.09). 
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Figure 3.24. Variation of normal force coefficient at 0.93 radius with HHC on, fixed-gain controller 
r = 0.0 (V = 135 knots, Cx/cr — 0.09). 
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Figure 3.25. Vibratory hub shear responses with fixed-gain HHC controller, r = 0.0 (V = 135 knots 
C T /a = 0.135). 
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Figure 3.26. Variation of normal force coefficient at 0.93 radius with HHC off (V — 135 knots, 
Cx/cr — 0.135). 
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Figure 3.27. Variation of normal force coefficient at 0.93 radius with HHC on, fixed-gain controller 
r = 0.0 ( V = 135 knots, C T /a = 0.135). 




Figure 3.28. Vibratory hub shear responses with adaptive-gain HHC controller, r = 0.3 ( V = 135 knots, 
C t /(7 = 0.135). 



Figure 3.29. Variation of normal force coefficient at 0.93 radius with HHC on, adaptive-gain controller 
r — 0.3 ( V = 135 knots, C^/cr = 0.135). 
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Figure 3.30. Stall area on rotor disk with HHC off (V = 135 knots, C T /a = 0.12). 








Figure 3.32. Variation of 3 P vertical hub shear with forward speed for different rotor systems 
(Ct/<j = 0.08). 



Figure 3.33. Variation of 2 P rotating inplane hub shear with forward speed for different rotor systems 
(C T /a = 0.08). J 
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Figure 3.34. Variation of AP rotating inplane hub shear with forward speed for different rotor systems 
(C T /a = 0.08). 7 



Figure 3.35. Variation of actuator power index with forward speed for different rotor systems 
(C T /cr = 0.08). 
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Figure 3.36. Effects of blade torsion stiffness on the uncontrolled vibratory hub shears (V = 40 knots 
Cj'/a — 0.08). ' 1 ' 5 



Figure 3.37. Effects of blade torsion stiffness on the uncontrolled vibratory hub moments 
(V = 40 knots, C T /cr = 0.08). 
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Figure 3.38. Effects of blade torsion stiffness on the HHC actuator power (V = 40 knots, Cj' / a = 0.08). 



Figure 3.39. Effects of blade flap stiffness on the uncontrolled vibratory hub shears (V = 40 knots 
C T /a = 0.08). 
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Figure 3.40. Effects of blade flap stiffness on the uncontrolled vibratory hub moments (V = 40 knots 
Cyycr = 0.08). 



Figure 3.41. Effects of blade flap stiffness on the HHC actuator power (V = 40 knots, Cj'/cr — 0.08). 
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Figure 3.42. Effects of blade lag stiffness on the uncontrolled vibratory hub shears (V 
C T /<7 = 0.08). 



Figure 3.43. Effects of blade lag stiffness on the uncontrolled vibratory hub moments (V 
C T /a = 0.08). 
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Figure 3.44. Effects of blade lag stiffness on the HHC actuator power (V — 40 knots, Cj'/tj = 0.08) 



Figure 3.45. Effects of blade linear pretwist on the uncontrolled vibratory hub shears (V = 40 knots 
C T /cr = 0.08). 
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Figure 3.46. Effects of blade linear pretwist on the uncontrolled vibratory hub moments (V = 40 knots 
C T /cr = 0.08). 
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Figure 3.47. Effects of blade linear pretwist on the HHC actuator power (V = 40 knots, C T /a = 0.08). 


113 






Figure 3 48. Effects of blade center-of-mass offset from elastic axis on the uncontrolled vibratory hub 
shears ( V = 40 knots, C T /cr = 0.08). 



Figure 3.49. Effects of blade center-of-mass offset from elastic axis on the uncontrolled vibratory hub 
moments (V = 40 knots, C T /a = 0.08). 
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Figure 3.50. Effects of blade center-of-mass offset from elastic axis on the HHC actuator power 
(V = 40 knots, C T / a = 0.08). 



Figure 3.51. Effects of blade elastic axis offset from quarter chord on the uncontrolled vibratory hub 
shears (V = 40 knots, C T /a = 0.08). 
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Figure 3.52. Effects of blade elastic axis offset from quarter chord on the uncontrolled vibratory hub 
moments (V = 40 knots, C T /a = 0.08). 
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Figure 3.53. Effects of blade elastic axis offset from quarter chord on the HHC actuator 
(V = 40 knots, CWcr = 0.08). 
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Figure 3.54. Effects of rotor thrust on the uncontrolled vibratory hub shears (V = 40 knots 
C T /a = 0.08). 
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Figure 3.55. Effects of rotor thrust on the uncontrolled vibratory hub moments ( V = 40 knots 
C T /a = 0.08). 
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Figure 3.56. Effects of rotor thrust on the HHC actuator power (V ~ 40 knots, Cx/cr = 0.08) 
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Figure 3.57. Effects of blade torsion stiffness on the uncontrolled vibratory hub shears ( V = 135 knots, 
C T /o- = 0 . 12 ). 
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BLADE TORSION FREQUENCY (PER REV) 

Figure 3.58. Effects of blade torsion stiffness on the HHC actuator power (V = 135 knots 
C T /cr = 0.12). 



FIRST FLAP BENDING FREQUENCY (PER REV) 

Figure 3.59. Effects of blade flap stiffness on the uncontrolled vibratory hub shears (V - 135 knots 
C T /a = 0.12). 
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FIRST FLAP BENDING FREQUENCY (PER REV) 

Figure 3.60. Effects of blade flap stiffness on the HHC actuator power (V = 135 knots, Cj>/a = 0.12). 



Figure 3.61. Effects of blade lag stiffness on the uncontrolled vibratory hub shears (V = 135 knots 
C T /a = 0.12). 
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FIRST LAG BENDING FREQUENCY (PER REV) 

Figure 3.62. Effects of blade lag stiffness on the HHC actuator power (V = 135 knots, Cj'j a — 0.12). 



Figure 3.63. Effects of blade center-of-mass offset from elastic axis on the uncontrolled vibratory hub 
shears {V = 135 knots, C T /a = 0.12). 
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Figure 3.64. Effects of blade center-of-mass offset from elastic axis on the HHC actuator power 
(V = 135 knots, C T /cr = 0.12). 
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Figure 3.65. Effects of blade elastic axis offset from quarter chord on the uncontrolled vibratory hub 
shears (V = 135 knots, C T /a = 0.12). 
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Figure 3.66. Effects of blade elastic axis offset from quarter chord on the HHC actuator power 
(V — 135 knots, Cx/c = 0.12). 
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